hankel API Summary¶
General quadrature method for Hankel transformations.
Based on the algorithm provided in H. Ogata, A Numerical Integration Formula Based on the Bessel Functions, Publications of the Research Institute for Mathematical Sciences, vol. 41, no. 4, pp. 949-970, 2005.
Classes¶
HankelTransform ([nu, N, h]) |
The basis of the Hankel Transformation algorithm by Ogata 2005. |
SymmetricFourierTransform ([ndim, a, b, N, h]) |
Determine the Fourier Transform of a radially symmetric function in arbitrary dimensions. |
Functions¶
get_h (f, nu[, K, cls, hstart, hdecrement, …]) |
Determine the largest value of h which gives a converged solution. |
-
class
hankel.
HankelTransform
(nu=0, N=None, h=0.05)[source]¶ The basis of the Hankel Transformation algorithm by Ogata 2005.
This algorithm is used to solve the equation \(\int_0^\infty f(x) J_\nu(x) dx\) where \(J_\nu(x)\) is a Bessel function of the first kind of order \(\nu\), and \(f(x)\) is an arbitrary (slowly-decaying) function.
The algorithm is presented in H. Ogata, A Numerical Integration Formula Based on the Bessel Functions, Publications of the Research Institute for Mathematical Sciences, vol. 41, no. 4, pp. 949-970, 2005.
This class provides a method for directly performing this integration, and also for doing a Hankel Transform.
Parameters: - nu : scalar, optional
The order of the bessel function (of the first kind) J_nu(x)
- N : int, optional, default =
pi/h
The number of nodes in the calculation. Generally this must increase for a smaller value of the step-size h. Default value is based on where the series will truncate according to the double-exponential convergence to the roots of the Bessel function.
- h : float, optional
The step-size of the integration.
-
classmethod
G
(f, h, k=None, *args, **kwargs)[source]¶ The absolute value of the non-oscillatory of the summed series’ last term, up to a scaling constant.
This can be used to get the sign of the slope of G with h.
Parameters: - f : callable
The function to integrate/transform
- h : float
The resolution parameter of the hankel integration
- k : float or array-like, optional
The scale at which to evaluate the transform. If None, assume an integral.
Returns: - The value of G.
-
classmethod
deltaG
(f, h, *args, **kwargs)[source]¶ The slope (up to a constant) of the last term of the series with h
-
integrate
(f, ret_err=True, ret_cumsum=False)[source]¶ Do the Hankel-type integral of the function f.
This is not the Hankel transform, but rather the simplified integral, \(\int_0^\infty f(x) J_\nu(x) dx\) , equivalent to the transform of \(f(r)/r\) at k=1.
Parameters: - f : callable
A function of one variable, representing \(f(x)\)
- ret_err : boolean, optional, default = True
Whether to return the estimated error
- ret_cumsum : boolean, optional, default = False
Whether to return the cumulative sum
-
transform
(f, k=1, ret_err=True, ret_cumsum=False, inverse=False)[source]¶ Do the Hankel-transform of the function f.
Parameters: - f : callable
A function of one variable, representing \(f(x)\)
- ret_err : boolean, optional, default = True
Whether to return the estimated error
- ret_cumsum : boolean, optional, default = False
Whether to return the cumulative sum
Returns: - ret : array-like
The Hankel-transform of f(x) at the provided k. If k is scalar, then this will be scalar.
- err : array-like
The estimated error of the approximate integral, at every k. It is merely the last term in the sum. Only returned if ret_err=True.
- cumsum : array-like
The total cumulative sum, for which the last term is itself the transform. One can use this to check whether the integral is converging. Only returned if ret_cumsum=True
Notes
The Hankel transform is defined as
\[F(k) = \int_0^\infty r f(r) J_\nu(kr) dr.\]The inverse transform is identical (swapping k and r of course).
-
xrange
(k=1)[source]¶ Tuple giving (min,max) x value evaluated by f(x).
Parameters: - k : array-like, optional
Scales for the transformation. Leave as 1 for an integral.
See also
See
- meth:xrange_approx for an approximate version of this method
which
-
classmethod
xrange_approx
(h, nu, k=1)[source]¶ Tuple giving approximate (min,max) x value evaluated by f(x/k).
Operates under the assumption that N = 3.2/h.
Parameters: - h : float
The resolution parameter of the Hankel integration
- nu : float
Order of the integration/transform
- k : array-like, optional
Scales for the transformation. Leave as 1 for an integral.
See also
xrange
- the actual x-range under a given choice of parameters.
-
class
hankel.
SymmetricFourierTransform
(ndim=2, a=1, b=1, N=200, h=0.05)[source]¶ Determine the Fourier Transform of a radially symmetric function in arbitrary dimensions.
Parameters: - ndim : int
Number of dimensions the transform is in.
- a, b : float, default 1
This pair of values defines the Fourier convention used (see Notes below for details)
- N : int, optional
The number of nodes in the calculation. Generally this must increase for a smaller value of the step-size h.
- h : float, optional
The step-size of the integration.
Notes
We allow for arbitrary Fourier convention, according to the scheme in http://mathworld.wolfram.com/FourierTransform.html. That is, we define the forward and inverse n-dimensional transforms respectively as
\[F(k) = \sqrt{\frac{|b|}{(2\pi)^{1-a}}}^n \int f(r) e^{i b\mathbf{k}\cdot\mathbf{r}} d^n\mathbf{r}\]and
\[f(r) = \sqrt{\frac{|b|}{(2\pi)^{1+a}}}^n \int F(k) e^{-i b\mathbf{k}\cdot\mathbf{r}} d^n \mathbf{k}.\]By default, we set both a and b to 1, so that the forward transform has a normalisation of unity.
In this general sense, the forward and inverse Hankel transforms are respectively
\[F(k) = \sqrt{\frac{|b|}{(2\pi)^{1-a}}}^n \frac{(2\pi)^{n/2}}{(bk)^{n/2-1}} \int_0^\infty r^{n/2-1} f(r) J_{n/2-1}(bkr) r dr\]and
\[f(r) = \sqrt{\frac{|b|}{(2\pi)^{1+a}}}^n \frac{(2\pi)^{n/2}}{(br)^{n/2-1}} \int_0^\infty k^{n/2-1} f(k) J_{n/2-1}(bkr) k dk.\]-
classmethod
G
(f, h, k=None, ndim=2)[source]¶ The absolute value of the non-oscillatory part of the summed series’ last term, up to a scaling constant.
This can be used to get the sign of the slope of G with h.
Parameters: - f : callable
The function to integrate/transform
- h : float
The resolution parameter of the hankel integration
- k : float or array-like, optional
The scale at which to evaluate the transform. If None, assume an integral.
- ndim : float
The number of dimensions of the transform
Returns: - The value of G.
-
transform
(f, k, *args, **kwargs)[source]¶ Do the n-symmetric Fourier transform of the function f.
Parameters and returns are precisely the same as
HankelTransform.transform()
.Notes
The n-symmetric fourier transform is defined in terms of the Hankel transform as
\[F(k) = \frac{(2\pi)^{n/2}}{k^{n/2-1}} \int_0^\infty r^{n/2-1} f(r) J_{n/2-1}(kr)r dr.\]The inverse transform has an inverse normalisation.
-
classmethod
xrange_approx
(h, ndim, k=1)[source]¶ Tuple giving approximate (min,max) x value evaluated by f(x/k).
Operates under the assumption that N = pi/h.
Parameters: - h : float
The resolution parameter of the Hankel integration
- ndim : float
Number of dimensions of the transform.
- k : array-like, optional
Scales for the transformation. Leave as 1 for an integral.
See also
xrange
- the actual x-range under a given choice of parameters.
-
hankel.
get_h
(f, nu, K=None, cls=<class 'hankel.hankel.HankelTransform'>, hstart=0.05, hdecrement=2, atol=0.001, rtol=0.001, maxiter=15, inverse=False)[source]¶ Determine the largest value of h which gives a converged solution.
Parameters: - f : callable
The function to be integrated/transformed.
- nu : float
Either the order of the transformation, or the number of dimensions (if cls is a
SymmetricFourierTransform
)- K : float or array-like, optional
The scale(s) of the transformation. If None, assumes an integration over f(x)J_nu(x) is desired. It is recommended to use a down-sampled K for this routine for efficiency. Often a min/max is enough.
- cls :
HankelTransform
subclass, optional Either
HankelTransform
or a subclass, specifying the type of transformation to do on f.- hstart : float, optional
The starting value of h.
- hdecrement : float, optional
How much to divide h by on each iteration.
- atol, rtol : float, optional
The tolerance parameters, passed to np.isclose, defining the stopping condition.
- maxiter : int, optional
Maximum number of iterations to perform.
- inverse : bool, optional
Whether to treat as an inverse transformation.
Returns: - h : float
The h value at which the solution converges.
- res : scalar or tuple
The value of the integral/transformation using the returned h – if a transformation, returns results at K.
- N : int
The number of nodes necessary in the final calculation. While each iteration uses N=pi/h, the returned N checks whether nodes are numerically zero above some threshold.
Notes
This function is not completely general. The function f is assumed to be reasonably smooth and non-oscillatory.
The idea is to use successively smaller values of h, with N=pi/h on each iteration, until the result betweeniterations becomes stable.