Python version - Methods

class glow.time_domain.It_SingleContour(Lens, y, p_prec={})

Bases: ItGeneral

Computation for the single contour regime (only one image).

Additional information: theory, default parameters.

(Only new parameters and attributes are documented. See ItGeneral for the internal information of the parent class)

Parameters:
p_precdict, optional

Precision parameters. New keys:

  • soft (float) – Softening factor.

  • method (str) – ODE integration method. Any kind recognized by Scipy is allowed, e.g. 'RK45'.

  • rtol (float) – Relative tolerance for the integration.

  • atol (float) – Absolute tolerance for the integration.

Methods

compute_novec(tau)

Compute \(I(\tau)\).

default_params()

(subclass defined) Initialize the default parameters.

deriv_phi(R, phi[, xc1, xc2])

Compute first derivatives of the Fermat potential wrt polar coordinates.

find_R_of_tau(tau)

Find the radius that corresponds to a given \(\tau\).

find_all_images()

Find the minimum of the Fermat potential assuming that it is the only critical point.

get_contour(tau[, N])

Contour corresponding to a given \(\tau\) in Cartesian coordinates.

integrate_contour(tau[, dense_output])

Integration of the contour corresponding to a given \(\tau\).

sys(t, u)

ODE system defining the contours (and \(I(t)\)).

x1x2(R, phi[, xc1, xc2])

Convert polar to Cartesian coordinates.

Warning

Compared to It_SingleContour_C, this class presents several drawbacks:

  • It does not check that we are indeed in the single contour regime. Proceed with caution…

  • The 'robust' option (i.e. parametric integration) is not implemented.

  • It is not parallelized and way slower than the C implementation.

default_params()

(subclass defined) Initialize the default parameters.

x1x2(R, phi, xc1=0, xc2=0)

Convert polar to Cartesian coordinates.

deriv_phi(R, phi, xc1=0, xc2=0)

Compute first derivatives of the Fermat potential wrt polar coordinates.

find_all_images()

Find the minimum of the Fermat potential assuming that it is the only critical point.

find_R_of_tau(tau)

Find the radius that corresponds to a given \(\tau\).

sys(t, u)

ODE system defining the contours (and \(I(t)\)).

integrate_contour(tau, dense_output=False)

Integration of the contour corresponding to a given \(\tau\).

Parameters:
taufloat

Relative time delay \(\tau\).

dense_outputbool

Include interpolation functions in the output.

Returns:
soloutput of Scipy’s solve_ivp

Solution of the differential equation.

compute_novec(tau)

Compute \(I(\tau)\). Non-vectorized version.

It includes the step function, \(I(\tau<0)=0\).

Parameters:
taufloat

Relative time delay \(\tau\).

Returns:
Ifloat

\(I(\tau)\)

get_contour(tau, N=1000)

Contour corresponding to a given \(\tau\) in Cartesian coordinates.

Parameters:
taufloat

Relative time delay \(\tau\).

Nint

Number of points in the contour.

Returns:
x1s, x2sarray

Contour in Cartesian coordinates.

class glow.time_domain.It_SingleIntegral(Lens, y, p_prec={})

Bases: ItGeneral

Computation for axisymmetric lenses, solving a single radial integral.

Additional information: theory, default parameters.

(Only new parameters and attributes are documented. See ItGeneral for the internal information of the parent class)

Parameters:
p_precdict, optional

Precision parameters. New keys:

  • method (str) – Function used for the radial integration. Options:

    • 'quad': Scipy’s quad, using QUADPACK.

    • 'quadrature': Scipy’s quadrature, fixed-tolerance Gaussian quadrature.

Methods

check_input()

(subclass defined, optional) Check the input of the implementation.

check_singcusp()

Check for the presence of singularities and cusps at the center.

classify_image(sol)

Classify an image.

compute(tau)

Computation of \(I(\tau)\).

default_params()

(subclass defined) Initialize the default parameters.

find_all_images()

Find all the critical points.

find_bracket(t)

Find all the different brackets of integration (i.e. all the different contours).

find_image_bracket_1d(x1_lo, x1_hi)

Find a root of \(\phi'(\boldsymbol{x})=0\), using a bracketing algorithm.

find_image_newton_1d(x1_guess)

Find a root of \(\phi'(\boldsymbol{x})=0\), using Newton's method.

find_moving_bracket_1d(xguess_lo, xguess_hi, t)

Find a root of \(\phi(x_1, 0)=t\), using a moving bracket method.

find_root_bracket_1d(x1_lo, x1_hi, t)

Find a root of \(\phi(x_1, 0)=t\), using a bracketing algorithm.

find_root_newton_1d(x1_guess, t)

Find a root of \(\phi(x_1, 0)=t\), using Newton's method.

get_contour(tau[, N])

Compute the contours in Cartesian coordinates.

integrate_quad(t)

Computation of \(I(t)\) with the 'quad' method.

integrate_quadrature(t)

Computation of \(I(t)\) with the 'quadrature' method.

Warning

The C version It_SingleIntegral_C should always be preferred. It uses the same algorithm, but the implementation is faster and more robust.

check_input()

(subclass defined, optional) Check the input of the implementation.

default_params()

(subclass defined) Initialize the default parameters.

find_image_newton_1d(x1_guess)

Find a root of \(\phi'(\boldsymbol{x})=0\), using Newton’s method.

Parameters:
x1_guessfloat

Initial guess.

Returns:
imagedict

Dictionary with the format explained in ItGeneral.

find_image_bracket_1d(x1_lo, x1_hi)

Find a root of \(\phi'(\boldsymbol{x})=0\), using a bracketing algorithm.

Parameters:
x1_lo, x1_hifloat

Bracket containing at least one root (function must change sign).

Returns:
imagedict

Dictionary with the format explained in ItGeneral.

find_root_newton_1d(x1_guess, t)

Find a root of \(\phi(x_1, 0)=t\), using Newton’s method.

Parameters:
x1_guessfloat

Initial guess.

tfloat

Time delay \(t\) (Note: related to the relative time delay as \(t=\tau+t_\text{min}\)).

Returns:
solScipy’s output of root_scalar

Solution.

find_root_bracket_1d(x1_lo, x1_hi, t)

Find a root of \(\phi(x_1, 0)=t\), using a bracketing algorithm.

Parameters:
x1_lo, x1_hifloat

Bracket containing at least one root (function must change sign).

tfloat

Time delay \(t\) (Note: related to the relative time delay as \(t=\tau+t_\text{min}\)).

Returns:
solScipy’s output of root_scalar

Solution.

find_moving_bracket_1d(xguess_lo, xguess_hi, t)

Find a root of \(\phi(x_1, 0)=t\), using a moving bracket method.

The algorithm starts with an initial bracket \([x_\text{lo}, x_\text{hi}]\). If \(f(x_\text{lo})f(x_\text{hi}) < 0\), there is a root in the bracket and a standard algorithm is used. Otherwise, the bracket is transformed to \([x_\text{lo}, x_\text{hi}]\to [x_\text{hi}, x_\text{hi}+\Delta x]\). The process is repeated, increasing \(\Delta x\) with each iteration, until the bracketing condition is met.

Parameters:
xguess_lo, xguess_hifloat

Initial bracket.

tfloat

Time delay \(t\) (Note: related to the relative time delay as \(t=\tau+t_\text{min}\)).

Returns:
solScipy’s output of root_scalar or None

Returns None if the bracketing condition is never met.

classify_image(sol)

Classify an image.

Parameters:
solScipy’s output of root_scalar

Solution to \(\phi'(\boldsymbol{x})=0\).

Returns:
imagedict

Dictionary with the image information. The format is explained in ItGeneral. If the solution has not converged, it returns an empty dictionary instead.

check_singcusp()

Check for the presence of singularities and cusps at the center.

Since \(\phi'(x=0)\) may diverge, we check if \(\partial_{x_1}\phi(-\delta x, 0)\phi(\delta x, 0) < 0\).

Returns:
singcuspdict

Image-like dictionary, with the format explained in ItGeneral.

find_all_images()

Find all the critical points.

Returns:
imageslist

Full list of images. The output is stored in p_crits, detailed in ItGeneral.

find_bracket(t)

Find all the different brackets of integration (i.e. all the different contours).

Parameters:
tfloat

Time delay \(t\).

Returns:
bracketslist

List of brackets with the limits of integration \([r_i, r_f]\).

integrate_quadrature(t)

Computation of \(I(t)\) with the 'quadrature' method.

Parameters:
tfloat or array

Time delay \(t\).

Returns:
Ifloat or array

\(I(t)\).

integrate_quad(t)

Computation of \(I(t)\) with the 'quad' method.

Parameters:
tfloat or array

Time delay \(t\).

Returns:
Ifloat or array

\(I(t)\).

compute(tau)

Computation of \(I(\tau)\).

Depending on the precision parameter method, it chooses between integrate_quadrature() or integrate_quad().

Parameters:
taufloat or array

Relative time delay \(\tau\).

Returns:
Ifloat or array

\(I(\tau)\).

get_contour(tau, N=1000)

Compute the contours in Cartesian coordinates.

Parameters:
taufloat

Relative time delay \(\tau\).

Nint

Number of points in the contour.

Returns:
contourlist

The length of this list is determined by the number of contours contributing to the given \(\tau\). Each element (i.e. each contour) is given in the format \([x_1, x_2]\).

class glow.time_domain.It_AnalyticSIS(y, p_prec={}, psi0=1)

Bases: ItGeneral

Analytic \(I(\tau)\) for the singular isothermal sphere.

Additional information: theory, default parameters.

Parameters:
psi0float

Normalization of the lensing potential \(\psi(x) = \psi_0 x\).

Methods

I_func(a, b, c, d[, soft])

Combination of complete elliptic integrals.

compute(tau)

Compute \(I(\tau)\).

default_params()

(subclass defined) Initialize the default parameters.

ellip_K(k)

Complete elliptic integral of the first kind, \(K(k)\).

ellip_Pi(n, k)

Complete elliptic integral of the third kind, \(\Pi(n, k)\).

find_images()

Compute (analytically) the properties of the images for the SIS.

default_params()

(subclass defined) Initialize the default parameters.

find_images()

Compute (analytically) the properties of the images for the SIS.

Returns:
imageslist

Full list of images. The output is stored in p_crits, detailed in ItGeneral.

ellip_K(k)

Complete elliptic integral of the first kind, \(K(k)\).

We follow the notation of [1].

Parameters:
kfloat or array
Returns:
Kfloat or array

References

[1]

I. S. Gradshteyn and I. M. Ryzhik, Table of Integrals, Series, and Products (Academic Press, 2007).

ellip_Pi(n, k)

Complete elliptic integral of the third kind, \(\Pi(n, k)\).

We follow the notation of [1].

Parameters:
n, kfloat or array
Returns:
Pifloat or array

References

[1]

I. S. Gradshteyn and I. M. Ryzhik, Table of Integrals, Series, and Products (Academic Press, 2007).

I_func(a, b, c, d, soft=0)

Combination of complete elliptic integrals.

\[I_\text{SIS}(a, b, c, d) = \frac{8(b-c)}{\sqrt{(a-c)(b-d)}} \left[ \Pi\left(\frac{a-b}{a-c}, r\right) + \frac{c\,K(r)}{(b-c)} \right]\]

where

\[r \equiv \sqrt{\frac{(a-b)(c-d)}{(a-c)(b-d)}}\]
Parameters:
a, b, c, dfloat or array
softfloat

Softening factor, to help avoid divergences in denominators.

Returns:
I_SISfloat or array
compute(tau)

Compute \(I(\tau)\).

Define the integral piece-wise in three regions:

  • Region 1: \((u > 1)\ \to\ (\tau > (\psi_0+y)^2/2)\)

  • Region 2: \((\sqrt{1-R^2} < u < 1)\ \to\ (2\psi_0 y < \tau < (\psi_0+y)^2/2)\)

    • 2A: \((R > 0)\ \to\ (\psi_0 > y)\)

    • 2B: \((R < 0)\ \to\ (\psi_0 < y)\)

  • Region 3: \((0 < u < \sqrt{1-R^2})\ \to\ (0 < \tau < 2\psi_0 y)\)

where

\[u \equiv \frac{\sqrt{2\tau}}{\psi_0 + y}\ ,\qquad R \equiv \frac{\psi_0-y}{\psi_0 + y}\]
Parameters:
taufloat or array

Relative time delay \(\tau\).

Returns:
Ifloat or array

\(I(\tau)\).

class glow.time_domain.It_NaiveAreaIntegral(Lens, y, p_prec={})

Bases: ItGeneral

Simple implementation of the binning/grid/area method for the computation of \(I(\tau)\).

The algorithm uses a uniform grid with polar coordinates centered at the lens position, \(\rho=r^2\) and \(\theta\). First, it stores all the evaluations of the Fermat potential into an array. This array is then sorted (thus we can directly compute the minimum time delay) and divided into \(\Delta t\)-boxes of unequal size. The size (i.e. temporal spacing) is chosen with the criterium that each box is computed using the same number of points.

This method is extraordinarily simple, provides the minimum time-delay and adaptative time sampling, but it cannot be scaled (since sorting becomes prohibitive). Also, with this method, some portion for large \(\tau\) must be thrown away (some of the contributions to the integral are always left out of the grid).

Additional information: theory, default parameters.

(Only new parameters and attributes are documented. See ItGeneral for the internal information of the parent class)

Parameters:
p_precdict, optional

Precision parameters. New keys:

  • rho_min, rho_max (float) – Minimum and maximum radial coordinate \(\rho=r^2\).

  • N_rho (int) – Number of points in the \(\rho\) axis.

  • N_theta (int) – Number of points in the \(\theta\) axis.

  • N_intervals (int) – Number of points in the evaluation of \(I(\tau)\).

Methods

compute()

Compute \(I(\tau)\).

default_params()

(subclass defined) Initialize the default parameters.

Warning

This method is only included due to its simplicity and lack of assumptions. It should never be used for applications that require high precision or speed. It is intended only for verification. The C version It_AreaIntegral_C provides a (slightly) better optimization, avoiding both the storage of all the grid evaluations and the time-sorting problem, using instead a fixed temporal grid.

default_params()

(subclass defined) Initialize the default parameters.

compute()

Compute \(I(\tau)\).

Returns:
t0float

Minimum time delay \(t_\text{min}\).

t_karray

Temporal grid \(\tau_i = t_i - t_\text{min}\).

It_karray

Result grid \(I_i=I(t_i)\).