# -*- coding: utf-8 -*-
r"""
Module Summary
--------------
Interfaces for the NAG Mark 30.3 `sum` Chapter.
``sum`` - Summation of Series
This module is concerned with the following tasks:
(a) calculating the **discrete Fourier transform** of a sequence of real or complex data values;
(#) calculating the **discrete convolution** or the **discrete correlation** of two sequences of real or complex data values using discrete Fourier transforms;
(#) calculating the **inverse Laplace transform** of a user-supplied function;
(#) calculating the **fast Gauss transform** approximation to the **discrete Gauss transform**;
(#) direct summation of orthogonal series;
(#) acceleration of convergence of a sequence of real values.
See Also
--------
``naginterfaces.library.examples.sum`` :
This subpackage contains examples for the ``sum`` module.
See also the :ref:`library_sum_ex` subsection.
Functionality Index
-------------------
Acceleration of convergence: :meth:`accelerate`
**Convolution or Correlation**
complex vectors: :meth:`convcorr_complex`
real vectors
time-saving: :meth:`convcorr_real`
**Discrete Fourier Transform**
multidimensional
complex sequence
complex storage: :meth:`fft_complex_multid`
real storage: :meth:`fft_complex_multid_sep`
multiple half - and quarter-wave transforms
Fourier cosine transforms
simple use: :meth:`fft_real_cosine_simple`
Fourier cosine transforms, simple use: :meth:`fft_cosine`
Fourier sine transforms
simple use: :meth:`fft_real_sine_simple`
Fourier sine transforms, simple use: :meth:`fft_sine`
quarter-wave cosine transforms
simple use: :meth:`fft_real_qtrcosine_simple`
quarter-wave cosine transforms, simple use: :meth:`fft_qtrcosine`
quarter-wave sine transforms
simple use: :meth:`fft_real_qtrsine_simple`
quarter-wave sine transforms, simple use: :meth:`fft_qtrsine`
one-dimensional
multiple transforms
complex sequence
complex storage, contiguous sequence elements: :meth:`fft_complex_1d_multi_col`
complex storage by rows: :meth:`fft_complex_1d_multi_row`
Hermitian/real sequence
complex storage, contiguous sequence elements: :meth:`fft_realherm_1d_multi_col`
complex storage, strided sequence elements: :meth:`fft_realherm_1d_multi_row`
multi-variable
complex sequence
complex storage: :meth:`fft_complex_multid_1d`
real storage: :meth:`fft_complex_multid_1d_sep`
single transforms
complex sequence
time-saving
complex storage: :meth:`fft_complex_1d`
real storage: :meth:`fft_complex_1d_sep`
Hermitian/real sequence
time-saving
complex storage: :meth:`fft_realherm_1d`
Hermitian sequence
time-saving
real storage: :meth:`fft_hermitian_1d_rfmt`
real sequence
time-saving
real storage: :meth:`fft_real_1d_rfmt`
three-dimensional
complex sequence
complex storage: :meth:`fft_complex_3d`
real storage: :meth:`fft_complex_3d_sep`
Hermitian/real sequence
complex-to-real: :meth:`fft_hermitian_3d`
real-to-complex: :meth:`fft_real_3d`
two-dimensional
complex sequence
complex storage: :meth:`fft_complex_2d`
Hermitian/real sequence
complex-to-real: :meth:`fft_hermitian_2d`
real-to-complex: :meth:`fft_real_2d`
Fast Gauss Transform: :meth:`fast_gauss`
**Inverse Laplace Transform**
Crump's method: :meth:`invlaplace_crump`
Weeks' method
compute coefficients of solution: :meth:`invlaplace_weeks`
evaluate solution: :meth:`invlaplace_weeks_eval`
Summation of Chebyshev series: :meth:`chebyshev`
For full information please refer to the NAG Library document
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06intro.html
"""
# NAG Copyright 2017-2024.
[docs]def accelerate(seqn, ncall, comm, nterms):
r"""
``accelerate`` accelerates the convergence of a given convergent sequence to its limit.
.. _c06ba-py2-py-doc:
For full information please refer to the NAG Library document for c06ba
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06baf.html
.. _c06ba-py2-py-parameters:
**Parameters**
**seqn** : float
The next term of the sequence to be considered.
**ncall** : int
On the first call :math:`\mathrm{ncall}` must be set to :math:`0`. Thereafter :math:`\mathrm{ncall}` must not be changed between calls.
**comm** : dict, communication object, modified in place
Communication structure.
`On initial entry`: need not be set.
**nterms** : int
The maximum number of terms in the sequence.
**Returns**
**ncall** : int
The number of terms in the sequence that have been considered.
**result** : float
The estimate of the limit of the sequence. For the first two calls, :math:`\mathrm{result} = \mathrm{seqn}`.
**abserr** : float
An estimate of the absolute error in :math:`\mathrm{result}`. For the first three calls, :math:`\mathrm{abserr}` is set to a large machine-dependent constant.
.. _c06ba-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{ncall} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{ncall} \geq 0`.
.. _c06ba-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``accelerate`` performs Shanks' transformation on a given sequence of real values by means of the Epsilon algorithm of Wynn (1956).
A (possibly unreliable) estimate of the absolute error is also given.
The function must be called repetitively, once for each new term in the sequence.
.. _c06ba-py2-py-references:
**References**
Shanks, D, 1955, `Nonlinear transformations of divergent and slowly convergent sequences`, J. Math. Phys. (34), 1--42
Wynn, P, 1956, `On a device for computing the` :math:`e_m\left(S_n\right)` `transformation`, Math. Tables Aids Comput. (10), 91--96
"""
raise NotImplementedError
[docs]def chebyshev(x, xmin, xmax, c, s):
r"""
``chebyshev`` evaluates a polynomial from its Chebyshev series representation at a set of points.
.. _c06dc-py2-py-doc:
For full information please refer to the NAG Library document for c06dc
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06dcf.html
.. _c06dc-py2-py-parameters:
**Parameters**
**x** : float, array-like, shape :math:`\left(\textit{lx}\right)`
:math:`x \in X`, the set of arguments of the series.
**xmin** : float
The lower end point of the interval :math:`\left[x_{\mathrm{min}}, x_{\mathrm{max}}\right]`.
**xmax** : float
The upper end point of the interval :math:`\left[x_{\mathrm{min}}, x_{\mathrm{max}}\right]`.
**c** : float, array-like, shape :math:`\left(n\right)`
:math:`\mathrm{c}[\textit{j}-1]` must contain the coefficient :math:`c_{\textit{j}}` of the Chebyshev series, for :math:`\textit{j} = 1,2,\ldots,n`.
**s** : int
Determines the series (see :ref:`Notes <c06dc-py2-py-notes>`).
:math:`\mathrm{s} = 1`
The series is general.
:math:`\mathrm{s} = 2`
The series is even.
:math:`\mathrm{s} = 3`
The series is odd.
**Returns**
**res** : float, ndarray, shape :math:`\left(\textit{lx}\right)`
The Chebyshev series evaluated at the set of points :math:`X`.
.. _c06dc-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\textit{lx} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\textit{lx} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{s} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{s} = 1`, :math:`2` or :math:`3`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{xmax} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{xmin} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{xmin} < \mathrm{xmax}`.
(`errno` :math:`5`)
On entry, element :math:`\mathrm{x}[\langle\mathit{\boldsymbol{value}}\rangle] = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{xmin} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{xmax} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{xmin}\leq \mathrm{x}[i]\leq \mathrm{xmax}`, for all :math:`i`.
.. _c06dc-py2-py-notes:
**Notes**
``chebyshev`` evaluates, at each point in a given set :math:`X`, the sum of a Chebyshev series of one of three forms according to the value of the parameter :math:`\mathrm{s}`:
.. rst-class:: nag-rules-none nag-align-left
+-----------------------+--------------------------------------------------------------+
|:math:`\mathrm{s} = 1`:|:math:`0.5c_1+\sum_{2}^{n}{c_j}T_{{j-1}}\left(\bar{x}\right)` |
+-----------------------+--------------------------------------------------------------+
|:math:`\mathrm{s} = 2`:|:math:`0.5c_1+\sum_{2}^{n}{c_j}T_{{2j-2}}\left(\bar{x}\right)`|
+-----------------------+--------------------------------------------------------------+
|:math:`\mathrm{s} = 3`:|:math:`\sum_{1}^{n}{c_j}T_{{2j-1}}\left(\bar{x}\right)` |
+-----------------------+--------------------------------------------------------------+
where :math:`\bar{x}` lies in the range :math:`{-1.0}\leq \bar{x}\leq 1.0`.
Here :math:`T_r\left(x\right)` is the Chebyshev polynomial of order :math:`r` in :math:`\bar{x}`, defined by :math:`\cos\left({ry}\right)` where :math:`\cos\left(y\right) = \bar{x}`.
It is assumed that the independent variable :math:`\bar{x}` in the interval :math:`\left[{-1.0}, {+1.0}\right]` was obtained from your original variable :math:`x \in X`, a set of real numbers in the interval :math:`\left[x_{\mathrm{min}}, x_{\mathrm{max}}\right]`, by the linear transformation
.. math::
\bar{x} = \frac{{2x-\left(x_{\mathrm{max}}+x_{\mathrm{min}}\right)}}{{x_{\mathrm{max}}-x_{\mathrm{min}}}}\text{.}
The method used is based upon a three-term recurrence relation; for details see Clenshaw (1962).
The coefficients :math:`c_j` are normally generated by other functions, for example they may be those returned by the interpolation function :meth:`interp.dim1_cheb <naginterfaces.library.interp.dim1_cheb>` (in vector :math:`\textit{a}`), by a least squares fitting function in submodule :mod:`~naginterfaces.library.fit`, or as the solution of a boundary value problem by :meth:`ode.bvp_coll_nth <naginterfaces.library.ode.bvp_coll_nth>`, :meth:`ode.bvp_coll_sys <naginterfaces.library.ode.bvp_coll_sys>` or :meth:`ode.bvp_ps_lin_solve <naginterfaces.library.ode.bvp_ps_lin_solve>`.
.. _c06dc-py2-py-references:
**References**
Clenshaw, C W, 1962, `Chebyshev Series for Mathematical Functions`, Mathematical tables, HMSO
"""
raise NotImplementedError
[docs]def fft_real_1d_rfmt(x):
r"""
``fft_real_1d_rfmt`` calculates the discrete Fourier transform of a sequence of :math:`n` real data values (using a work array for extra speed).
.. _c06fa-py2-py-doc:
For full information please refer to the NAG Library document for c06fa
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06faf.html
.. _c06fa-py2-py-parameters:
**Parameters**
**x** : float, array-like, shape :math:`\left(n\right)`
If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_real_1d_rfmt`` is called, :math:`\mathrm{x}[\textit{j}]` must contain :math:`x_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n\right)`
The discrete Fourier transform stored in Hermitian form. If the components of the transform :math:`\hat{z}_k` are written as :math:`a_k+ib_k`, and if :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_real_1d_rfmt`` is called, then for :math:`0\leq k\leq {n/2}`, :math:`a_k` is contained in :math:`\mathrm{x}[k-1]`, and for :math:`1\leq k\leq {\left(n-1\right)/2}`, :math:`b_k` is contained in :math:`\mathrm{x}[n-k]`. (See also `the C06 Introduction <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06intro.html#background12>`__.)
.. _c06fa-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`3`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n > 1`.
.. _c06fa-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given a sequence of :math:`n` real data values :math:`x_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_real_1d_rfmt`` calculates their discrete Fourier transform defined by
.. math::
\hat{z}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}x_j\times \mathrm{exp}\left({-i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{.}
(Note the scale factor of :math:`\frac{1}{\sqrt{n}}` in this definition.) The transformed values :math:`\hat{z}_k` are complex, but they form a Hermitian sequence (i.e., :math:`\hat{z}_{{n-k}}` is the complex conjugate of :math:`\hat{z}_k`), so they are completely determined by :math:`n` real numbers (see also `the C06 Introduction <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06intro.html>`__).
To compute the inverse discrete Fourier transform defined by
.. math::
\hat{w}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}x_j\times \mathrm{exp}\left({+i}\frac{{2\pi jk}}{n}\right)\text{,}
this function should be followed by forming the complex conjugates of the :math:`\hat{z}_k`; that is, :math:`x\left(\textit{k}\right) = -x\left(\textit{k}\right)`, for :math:`\textit{k} = n/2+2,\ldots,n`.
``fft_real_1d_rfmt`` uses the fast Fourier transform (FFT) algorithm (see Brigham (1974)).
.. _c06fa-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
"""
raise NotImplementedError
[docs]def fft_hermitian_1d_rfmt(x):
r"""
``fft_hermitian_1d_rfmt`` calculates the discrete Fourier transform of a Hermitian sequence of :math:`n` complex data values (using a work array for extra speed).
.. _c06fb-py2-py-doc:
For full information please refer to the NAG Library document for c06fb
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06fbf.html
.. _c06fb-py2-py-parameters:
**Parameters**
**x** : float, array-like, shape :math:`\left(n\right)`
The sequence to be transformed stored in Hermitian form. If the data values :math:`z_j` are written as :math:`x_j+iy_j`, and if :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_hermitian_1d_rfmt`` is called, then for :math:`0\leq j\leq {n/2}`, :math:`x_j` is contained in :math:`\mathrm{x}[j-1]`, and for :math:`1\leq j\leq {\left(n-1\right)/2}`, :math:`y_j` is contained in :math:`\mathrm{x}[n-j]`. (See also `the C06 Introduction <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06intro.html#background12>`__.)
**Returns**
**x** : float, ndarray, shape :math:`\left(n\right)`
The components of the discrete Fourier transform :math:`\hat{x}_k`. If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_hermitian_1d_rfmt`` is called, :math:`\hat{x}_{\textit{k}}` is stored in :math:`\mathrm{x}[\textit{k}]`, for :math:`\textit{k} = 0,1,\ldots,{n-1}`.
.. _c06fb-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`3`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n > 1`.
.. _c06fb-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given a Hermitian sequence of :math:`n` complex data values :math:`z_{\textit{j}}` (i.e., a sequence such that :math:`z_0` is real and :math:`z_{{n-\textit{j}}}` is the complex conjugate of :math:`z_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,{n-1}`), ``fft_hermitian_1d_rfmt`` calculates their discrete Fourier transform defined by
.. math::
\hat{x}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j\times \mathrm{exp}\left({-i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{.}
(Note the scale factor of :math:`\frac{1}{\sqrt{n}}` in this definition.) The transformed values :math:`\hat{x}_k` are purely real (see also `the C06 Introduction <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06intro.html>`__).
To compute the inverse discrete Fourier transform defined by
.. math::
\hat{y}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j\times \mathrm{exp}\left({+i}\frac{{2\pi jk}}{n}\right)\text{,}
this function should be preceded by forming the complex conjugates of the :math:`\hat{z}_k`; that is, :math:`x\left(\textit{k}\right) = -x\left(\textit{k}\right)`, for :math:`\textit{k} = n/2+2,\ldots,n`.
``fft_hermitian_1d_rfmt`` uses the fast Fourier transform (FFT) algorithm (see Brigham (1974)).
.. _c06fb-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
"""
raise NotImplementedError
[docs]def fft_complex_1d_sep(x, y):
r"""
``fft_complex_1d_sep`` calculates the discrete Fourier transform of a sequence of :math:`n` complex data values (using a work array for extra speed).
.. _c06fc-py2-py-doc:
For full information please refer to the NAG Library document for c06fc
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06fcf.html
.. _c06fc-py2-py-parameters:
**Parameters**
**x** : float, array-like, shape :math:`\left(n\right)`
If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_complex_1d_sep`` is called, :math:`\mathrm{x}[\textit{j}-1]` must contain :math:`x_{\textit{j}}`, the real part of :math:`z_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
**y** : float, array-like, shape :math:`\left(n\right)`
If :math:`\mathrm{y}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_complex_1d_sep`` is called, :math:`\mathrm{y}[\textit{j}-1]` must contain :math:`y_{\textit{j}}`, the imaginary part of :math:`z_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n\right)`
The real parts :math:`a_k` of the components of the discrete Fourier transform. If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_complex_1d_sep`` is called, for :math:`0\leq k\leq {n-1}`, :math:`a_k` is contained in :math:`\mathrm{x}[k-1]`.
**y** : float, ndarray, shape :math:`\left(n\right)`
The imaginary parts :math:`b_k` of the components of the discrete Fourier transform. If :math:`\mathrm{y}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_complex_1d_sep`` is called, then for :math:`0\leq k\leq {n-1}`, :math:`b_k` is contained in :math:`\mathrm{y}[k-1]`.
.. _c06fc-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`3`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n > 1`.
.. _c06fc-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given a sequence of :math:`n` complex data values :math:`z_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_complex_1d_sep`` calculates their discrete Fourier transform defined by
.. math::
\hat{z}_k = a_k+ib_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j\times \mathrm{exp}\left({-i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{.}
(Note the scale factor of :math:`\frac{1}{\sqrt{n}}` in this definition.)
To compute the inverse discrete Fourier transform defined by
.. math::
\hat{w}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j\times \mathrm{exp}\left({+i}\frac{{2\pi jk}}{n}\right)\text{,}
this function should be preceded and followed by the complex conjugation of the data values and the transform (by negating the imaginary parts stored in :math:`y`).
``fft_complex_1d_sep`` uses the fast Fourier transform (FFT) algorithm (see Brigham (1974)).
.. _c06fc-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
"""
raise NotImplementedError
[docs]def fft_complex_multid_1d_sep(l, nd, x, y):
r"""
``fft_complex_multid_1d_sep`` computes the discrete Fourier transform of one variable in a multivariate sequence of complex data values.
.. _c06ff-py2-py-doc:
For full information please refer to the NAG Library document for c06ff
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06fff.html
.. _c06ff-py2-py-parameters:
**Parameters**
**l** : int
:math:`l`, the index of the variable (or dimension) on which the discrete Fourier transform is to be performed.
**nd** : int, array-like, shape :math:`\left(\textit{ndim}\right)`
:math:`\mathrm{nd}[\textit{i}-1]` must contain :math:`n_{\textit{i}}` (the dimension of the :math:`\textit{i}`\ th variable), for :math:`\textit{i} = 1,2,\ldots,m`.
**x** : float, array-like, shape :math:`\left(n\right)`
:math:`\mathrm{x}[ 1 + j_1 + n_1 j_2 + n_1 n_2 j_3 + \cdots -1]` must contain the real part of the complex data value :math:`z_{{j_1j_2\ldots j_m}}`, for :math:`0\leq j_1\leq n_1-1,0\leq j_2\leq n_2-1,\ldots \text{}`; i.e., the values are stored in consecutive elements of the array according to the Fortran convention for storing multidimensional arrays.
**y** : float, array-like, shape :math:`\left(n\right)`
The imaginary parts of the complex data values, stored in the same way as the real parts in the array :math:`\mathrm{x}`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n\right)`
The real parts of the corresponding elements of the computed transform.
**y** : float, ndarray, shape :math:`\left(n\right)`
The imaginary parts of the corresponding elements of the computed transform.
.. _c06ff-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\textit{ndim} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\textit{ndim} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n = {\mathrm{nd}[0]\times \mathrm{nd}[1]\times \cdots \times \mathrm{nd}[\textit{ndim}-1]}`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{l} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\textit{ndim} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`1\leq \mathrm{l}\leq \textit{ndim}`.
(`errno` :math:`{10\times \mathrm{l}+3}`)
On entry, :math:`l = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{nd}[l-1] \geq 1`.
.. _c06ff-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``fft_complex_multid_1d_sep`` computes the discrete Fourier transform of one variable (the :math:`l`\ th say) in a multivariate sequence of complex data values :math:`z_{{j_1j_2\ldots j_m}}`, for :math:`\textit{j}_2 = 0,1,\ldots,{n_2-1}`, for :math:`\textit{j}_1 = 0,1,\ldots,{n_1-1}`, and so on.
Thus the individual dimensions are :math:`n_1,n_2,\ldots,n_m`, and the total number of data values is :math:`n = n_1\times n_2\times \cdots \times n_m`.
The function computes :math:`n/n_l` one-dimensional transforms defined by
.. math::
\hat{z}_{{j_1\ldots k_l\ldots j_m}} = \frac{1}{\sqrt{n_l}}\sum_{{j_l = 0}}^{{n_l-1}}z_{{j_1\ldots j_l\ldots j_m}}\times \mathrm{exp}\left(-\frac{{2\pi ij_lk_l}}{n_l}\right)\text{,}
where :math:`k_l = 0,1,\ldots,n_l-1`.
(Note the scale factor of :math:`\frac{1}{\sqrt{n_l}}` in this definition.)
To compute the inverse discrete Fourier transforms, defined with :math:`\mathrm{exp}\left(+\frac{{2\pi ij_lk_l}}{n_l}\right)` in the above formula instead of :math:`\mathrm{exp}\left(-\frac{{2\pi ij_lk_l}}{n_l}\right)`, this function should be preceded and followed by the complex conjugation of the data values and the transform (by negating the imaginary parts stored in :math:`y`).
The data values must be supplied in a pair of one-dimensional arrays (real and imaginary parts separately), in accordance with the Fortran convention for storing multidimensional data (i.e., with the first subscript :math:`j_1` varying most rapidly).
This function calls :meth:`fft_complex_1d_sep` to perform one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm in Brigham (1974).
.. _c06ff-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
"""
raise NotImplementedError
[docs]def fft_complex_multid_sep(nd, x, y):
r"""
``fft_complex_multid_sep`` computes the multidimensional discrete Fourier transform of a multivariate sequence of complex data values.
.. _c06fj-py2-py-doc:
For full information please refer to the NAG Library document for c06fj
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06fjf.html
.. _c06fj-py2-py-parameters:
**Parameters**
**nd** : int, array-like, shape :math:`\left(\textit{ndim}\right)`
:math:`\mathrm{nd}[\textit{i}-1]` must contain :math:`n_{\textit{i}}` (the dimension of the :math:`\textit{i}`\ th variable), for :math:`\textit{i} = 1,2,\ldots,m`.
**x** : float, array-like, shape :math:`\left(n\right)`
:math:`\mathrm{x}[ 1 + j_1 + n_1 j_2 + n_1 n_2 j_3 + \cdots -1]` must contain the real part of the complex data value :math:`z_{{j_1j_2\ldots j_m}}`, for :math:`0\leq j_1\leq n_1-1,0\leq j_2\leq n_2-1,\ldots \text{}`; i.e., the values are stored in consecutive elements of the array according to the Fortran convention for storing multidimensional arrays.
**y** : float, array-like, shape :math:`\left(n\right)`
The imaginary parts of the complex data values, stored in the same way as the real parts in the array :math:`\mathrm{x}`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n\right)`
The real parts of the corresponding elements of the computed transform.
**y** : float, ndarray, shape :math:`\left(n\right)`
The imaginary parts of the corresponding elements of the computed transform.
.. _c06fj-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\textit{ndim} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\textit{ndim} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n = {\mathrm{nd}[0]\times \mathrm{nd}[1]\times \cdots \times \mathrm{nd}[\textit{ndim}-1]}`.
(`errno` :math:`i \bmod 10 = 3`)
On entry, :math:`l = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{nd}[l-1] \geq 1`.
.. _c06fj-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``fft_complex_multid_sep`` computes the multidimensional discrete Fourier transform of a multidimensional sequence of complex data values :math:`z_{{j_1j_2\ldots j_m}}`, where :math:`j_1 = 0,1,\ldots,n_1-1\text{, }\quad j_2 = 0,1,\ldots,n_2-1`, and so on.
Thus the individual dimensions are :math:`n_1,n_2,\ldots,n_m`, and the total number of data values is :math:`n = n_1\times n_2\times \cdots \times n_m`.
The discrete Fourier transform is here defined (e.g., for :math:`m = 2`) by:
.. math::
\hat{z}_{{k_1,k_2}} = \frac{1}{\sqrt{n}}\sum_{{j_1 = 0}}^{{n_1-1}}\sum_{{j_2 = 0}}^{{n_2-1}}z_{{j_1j_2}}\times \mathrm{exp}\left({-2\pi i}\left(\frac{{j_1k_1}}{n_1}+\frac{{j_2k_2}}{n_2}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,n_1-1`, :math:`k_2 = 0,1,\ldots,n_2-1`.
The extension to higher dimensions is obvious. (Note the scale factor of :math:`\frac{1}{\sqrt{n}}` in this definition.)
To compute the inverse discrete Fourier transform, defined with :math:`\mathrm{exp}\left(+2\pi i\left(\ldots \right)\right)` in the above formula instead of :math:`\mathrm{exp}\left(-2\pi i\left(\ldots \right)\right)`, this function should be preceded and followed by the complex conjugation of the data values and the transform (by negating the imaginary parts stored in :math:`y`).
The data values must be supplied in a pair of one-dimensional arrays (real and imaginary parts separately), in accordance with the Fortran convention for storing multidimensional data (i.e., with the first subscript :math:`j_1` varying most rapidly).
This function calls :meth:`fft_complex_1d_sep` to perform one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm in Brigham (1974).
.. _c06fj-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
"""
raise NotImplementedError
[docs]def convcorr_real(job, x, y):
r"""
``convcorr_real`` calculates the circular convolution or correlation of two real vectors of period :math:`n` (using a work array for extra speed).
.. _c06fk-py2-py-doc:
For full information please refer to the NAG Library document for c06fk
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06fkf.html
.. _c06fk-py2-py-parameters:
**Parameters**
**job** : int
The computation to be performed.
:math:`\mathrm{job} = 1`
:math:`z_k = {\sum_{{j = 0}}^{{n-1}}x_jy_{{k-j}}}` (convolution);
:math:`\mathrm{job} = 2`
:math:`w_k = {\sum_{{j = 0}}^{{n-1}}x_jy_{{k+j}}}` (correlation).
**x** : float, array-like, shape :math:`\left(n\right)`
The elements of one period of the vector :math:`x`. If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``convcorr_real`` is called, :math:`\mathrm{x}[\textit{j}]` must contain :math:`x_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
**y** : float, array-like, shape :math:`\left(n\right)`
The elements of one period of the vector :math:`y`. If :math:`\mathrm{y}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``convcorr_real`` is called, :math:`\mathrm{y}[\textit{j}]` must contain :math:`y_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n\right)`
The corresponding elements of the discrete convolution or correlation.
**y** : float, ndarray, shape :math:`\left(n\right)`
The discrete Fourier transform of the convolution or correlation returned in the array :math:`\mathrm{x}`; the transform is stored in Hermitian form; if the components of the transform :math:`z_k` are written as :math:`a_k+ib_k`, then for :math:`0\leq k\leq n/2`, :math:`a_k` is contained in :math:`\mathrm{y}[k]`, and for :math:`1\leq k\leq n/2-1`, :math:`b_k` is contained in :math:`\mathrm{y}[n-k]`. (See also `the C06 Introduction <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06intro.html#background12>`__.)
.. _c06fk-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`3`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n > 1`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{job} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{job} = 1` or :math:`2`.
.. _c06fk-py2-py-notes:
**Notes**
``convcorr_real`` computes:
if :math:`\mathrm{job} = 1`, the discrete **convolution** of :math:`x` and :math:`y`, defined by
.. math::
z_k = \sum_{{j = 0}}^{{n-1}}x_jy_{{k-j}} = \sum_{{j = 0}}^{{n-1}}x_{{k-j}}y_j\text{;}
if :math:`\mathrm{job} = 2`, the discrete **correlation** of :math:`x` and :math:`y` defined by
.. math::
w_k = \sum_{{j = 0}}^{{n-1}}x_jy_{{k+j}}\text{.}
Here :math:`x` and :math:`y` are real vectors, assumed to be periodic, with period :math:`n`, i.e., :math:`x_j = x_{{j\pm n}} = x_{{j\pm 2n}} = \cdots \text{}`; :math:`z` and :math:`w` are then also periodic with period :math:`n`.
**Note:** this usage of the terms 'convolution' and 'correlation' is taken from Brigham (1974).
The term 'convolution' is sometimes used to denote both these computations.
If :math:`\hat{x}`, :math:`\hat{y}`, :math:`\hat{z}` and :math:`\hat{w}` are the discrete Fourier transforms of these sequences, i.e.,
.. math::
\hat{x}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}x_j\times \mathrm{exp}\left({-i}\frac{{2\pi jk}}{n}\right)\text{, etc.,}
then :math:`\hat{z}_k = \sqrt{n}.\hat{x}_k\hat{y}_k` and :math:`\hat{w}_k = \sqrt{n}.\bar{\hat{x}}_k\hat{y}_k` (the bar denoting complex conjugate).
This function calls the same auxiliary functions as :meth:`fft_realherm_1d` to compute discrete Fourier transforms.
.. _c06fk-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
See Also
--------
:meth:`naginterfaces.library.examples.sum.convcorr_real_ex.main`
"""
raise NotImplementedError
[docs]def fft_complex_3d_sep(n1, n2, n3, x, y):
r"""
``fft_complex_3d_sep`` computes the three-dimensional discrete Fourier transform of a trivariate sequence of complex data values.
This function is designed to be particularly efficient on vector processors.
.. _c06fx-py2-py-doc:
For full information please refer to the NAG Library document for c06fx
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06fxf.html
.. _c06fx-py2-py-parameters:
**Parameters**
**n1** : int
:math:`n_1`, the first dimension of the transform.
**n2** : int
:math:`n_2`, the second dimension of the transform.
**n3** : int
:math:`n_3`, the third dimension of the transform.
**x** : float, array-like, shape :math:`\left(\mathrm{n1}\times \mathrm{n2}\times \mathrm{n3}\right)`
The real and imaginary parts of the complex data values must be stored in arrays :math:`\mathrm{x}` and :math:`\mathrm{y}` respectively. If :math:`\mathrm{x}` and :math:`\mathrm{y}` are regarded as three-dimensional arrays of dimension :math:`\left({0:\mathrm{n1}-1}, {0:\mathrm{n2}-1}, {0:\mathrm{n3}-1}\right)`, :math:`\mathrm{x}[j_1-1,j_2-1,j_3-1]` and :math:`\mathrm{y}[j_1-1,j_2-1,j_3-1]` must contain the real and imaginary parts of :math:`z_{{j_1j_2j_3}}`.
**y** : float, array-like, shape :math:`\left(\mathrm{n1}\times \mathrm{n2}\times \mathrm{n3}\right)`
The real and imaginary parts of the complex data values must be stored in arrays :math:`\mathrm{x}` and :math:`\mathrm{y}` respectively. If :math:`\mathrm{x}` and :math:`\mathrm{y}` are regarded as three-dimensional arrays of dimension :math:`\left({0:\mathrm{n1}-1}, {0:\mathrm{n2}-1}, {0:\mathrm{n3}-1}\right)`, :math:`\mathrm{x}[j_1-1,j_2-1,j_3-1]` and :math:`\mathrm{y}[j_1-1,j_2-1,j_3-1]` must contain the real and imaginary parts of :math:`z_{{j_1j_2j_3}}`.
**Returns**
**x** : float, ndarray, shape :math:`\left(\mathrm{n1}\times \mathrm{n2}\times \mathrm{n3}\right)`
The real and imaginary parts respectively of the corresponding elements of the computed transform.
**y** : float, ndarray, shape :math:`\left(\mathrm{n1}\times \mathrm{n2}\times \mathrm{n3}\right)`
The real and imaginary parts respectively of the corresponding elements of the computed transform.
.. _c06fx-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{n1} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n1} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n2} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n2} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{n3} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n3} \geq 1`.
.. _c06fx-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``fft_complex_3d_sep`` computes the three-dimensional discrete Fourier transform of a trivariate sequence of complex data values :math:`z_{{j_1j_2j_3}}`, for :math:`\textit{j}_3 = 0,1,\ldots,{n_3-1}`, for :math:`\textit{j}_2 = 0,1,\ldots,{n_2-1}`, for :math:`\textit{j}_1 = 0,1,\ldots,{n_1-1}`.
The discrete Fourier transform is here defined by
.. math::
\hat{z}_{{k_1k_2k_3}} = \frac{1}{\sqrt{n_1n_2n_3}}\sum_{{j_1 = 0}}^{{n_1-1}}\sum_{{j_2 = 0}}^{{n_2-1}}\sum_{{j_3 = 0}}^{{n_3-1}}z_{{j_1j_2j_3}}\times \mathrm{exp}\left({-2\pi i}\left(\frac{{j_1k_1}}{n_1}+\frac{{j_2k_2}}{n_2}+\frac{{j_3k_3}}{n_3}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,{n_1-1}`, :math:`k_2 = 0,1,\ldots,{n_2-1}`, :math:`k_3 = 0,1,\ldots,{n_3-1}`.
(Note the scale factor of :math:`\frac{1}{\sqrt{n_1n_2n_3}}` in this definition.)
To compute the inverse discrete Fourier transform, defined with :math:`\mathrm{exp}\left({+2\pi i}\left(\ldots \right)\right)` in the above formula instead of :math:`\mathrm{exp}\left({-2\pi i}\left(\ldots \right)\right)`, this function should be preceded and followed by forming the complex conjugates of the data values and the transform.
This function performs, for each dimension, multiple one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm (see Brigham (1974)).
It is designed to be particularly efficient on vector processors.
.. _c06fx-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def invlaplace_crump(f, t, relerr, alphab, tfac=0.8, mxterm=100, data=None):
r"""
``invlaplace_crump`` estimates values of the inverse Laplace transform of a given function using a Fourier series approximation. A bound on the exponential order of the inverse is required.
.. _c06la-py2-py-doc:
For full information please refer to the NAG Library document for c06la
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06laf.html
.. _c06la-py2-py-parameters:
**Parameters**
**f** : callable fp = f(p, data=None)
:math:`\mathrm{f}` must evaluate the function :math:`F\left(p\right)` for a given value of :math:`p`.
**Parameters**
**p** : complex
The argument :math:`p`.
**data** : arbitrary, optional, modifiable in place
User-communication data for callback functions.
**Returns**
**fp** : complex
The value :math:`F\left(p\right)`.
**t** : float, array-like, shape :math:`\left(n\right)`
Each :math:`\mathrm{t}[\textit{j}-1]` must specify a point at which the inverse Laplace transform is required, for :math:`\textit{j} = 1,2,\ldots,n`.
**relerr** : float
The required relative error in the values of the inverse Laplace transform. If the absolute value of the inverse is less than :math:`\mathrm{relerr}`, absolute accuracy is used instead. :math:`\mathrm{relerr}` must be in the range :math:`0.0\leq \mathrm{relerr} < 1.0`. If :math:`\mathrm{relerr}` is set too small or to :math:`0.0`, the function uses a value sufficiently larger than machine precision.
**alphab** : float
:math:`\alpha_b`, an upper bound for :math:`\alpha` (see :ref:`Notes <c06la-py2-py-notes>`). Usually, :math:`\alpha_b` should be specified equal to, or slightly larger than, the value of :math:`\alpha`. If :math:`\alpha_b < \alpha` then the prescribed accuracy may not be achieved or completely incorrect results may be obtained. If :math:`\alpha_b` is too large ``invlaplace_crump`` will be inefficient and convergence may not be achieved.
**Note:** it is as important to specify :math:`\alpha_b` correctly as it is to specify the correct function for inversion.
**tfac** : float, optional
:math:`t_{\mathrm{fac}}`, a factor to be used in calculating the parameter :math:`\tau`. Larger values (e.g., :math:`5.0`) may be specified for difficult problems, but these may require very large values of :math:`\mathrm{mxterm}`.
**mxterm** : int, optional
The maximum number of (complex) terms to be used in the evaluation of the Fourier series.
**data** : arbitrary, optional
User-communication data for callback functions.
**Returns**
**valinv** : float, ndarray, shape :math:`\left(n\right)`
An estimate of the value of the inverse Laplace transform at :math:`t = \mathrm{t}[\textit{j}-1]`, for :math:`\textit{j} = 1,2,\ldots,n`.
**errest** : float, ndarray, shape :math:`\left(n\right)`
An estimate of the error in :math:`\mathrm{valinv}[j]`. This is usually an estimate of relative error but, if :math:`\mathrm{valinv}[j-1] < \mathrm{relerr}`, :math:`\mathrm{errest}[j]` estimates the absolute error. :math:`\mathrm{errest}[j-1]` is unreliable when :math:`\mathrm{valinv}[j-1]` is small but slightly greater than :math:`\mathrm{relerr}`.
**nterms** : int
The number of (complex) terms actually used.
**na** : int
The number of values of :math:`a` used by the function. See `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06laf.html#fcomments>`__.
**alow** : float
The smallest value of :math:`a` used in the algorithm. This may be used for checking the value of :math:`\mathrm{alphab}-` see `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06laf.html#fcomments>`__.
**ahigh** : float
The largest value of :math:`a` used in the algorithm. This may be used for checking the value of :math:`\mathrm{alphab}-` see `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06laf.html#fcomments>`__.
**nfeval** : int
The number of calls to :math:`\mathrm{f}` made by the function.
.. _c06la-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{tfac} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{tfac} > 0.5`.
(`errno` :math:`1`)
On entry, :math:`\mathrm{relerr} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`0.0\leq \mathrm{relerr} < 1.0`.
(`errno` :math:`1`)
On entry, :math:`\mathrm{mxterm} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{mxterm} \geq 1`.
(`errno` :math:`1`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`2`)
On entry, the elements of :math:`\mathrm{t}` are not strictly increasing.
(`errno` :math:`2`)
On entry, :math:`\mathrm{t}[0] = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`0.0\leq \mathrm{t}[0] < \mathrm{t}[1] < \cdots < \mathrm{t}[n-1]`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{t}[n-1]` is too large: :math:`\mathrm{t}[n-1] = \langle\mathit{\boldsymbol{value}}\rangle`. If necessary, scale the problem as described in `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06laf.html#fcomments>`__.
(`errno` :math:`4`)
Required accuracy cannot be obtained. Try increasing :math:`\mathrm{tfac}`, :math:`\mathrm{alphab}`, or both. :math:`\mathrm{tfac} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{alphab} = \langle\mathit{\boldsymbol{value}}\rangle`.
(`errno` :math:`5`)
Convergence failure in epsilon algorithm: :math:`\mathrm{mxterm} = \langle\mathit{\boldsymbol{value}}\rangle`. Some values of :math:`\mathrm{valinv}[j-1]` may be calculated to the desired accuracy; this may be determined by examining the values of :math:`\mathrm{errest}[j-1]`. Try reducing the range of :math:`\mathrm{t}` or increasing :math:`\mathrm{mxterm}`. If :math:`\mathrm{errno}` = 5 still results, try reducing :math:`\mathrm{tfac}`.
**Warns**
**NagAlgorithmicWarning**
(`errno` :math:`6`)
All values of :math:`\mathrm{valinv}` have been calculated, but are not of the required accuracy; the values of :math:`\mathrm{errest}[j-1]` should be examined carefully. Try reducing the range of :math:`\mathrm{t}`, or increasing :math:`\mathrm{tfac}`, :math:`\mathrm{alphab}` or both.
.. _c06la-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given a function :math:`F\left(p\right)` defined for complex values of :math:`p`, ``invlaplace_crump`` estimates values of its inverse Laplace transform by Crump's method (see Crump (1976)). (For a definition of the Laplace transform and its inverse, see `the C06 Introduction <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06intro.html>`__.)
Crump's method applies the epsilon algorithm (see Wynn (1956)) to the summation in Durbin's Fourier series approximation (see Durbin (1974))
.. math::
f\left(t_j\right)\simeq \frac{e^{{at_j}}}{\tau }\left[\frac{1}{2}F\left(a\right)-\sum_{{k = 1}}^{\infty }\left\{\mathrm{Re}\left(F\left(a+\frac{{k\pi i}}{\tau }\right)\right)\cos\left(\frac{{k\pi t_j}}{\tau }\right)-\mathrm{Im}\left(F\left(a+\frac{{k\pi i}}{\tau }\right)\right)\sin\left(\frac{{k\pi t_j}}{\tau }\right)\right\}\right]\text{,}
for :math:`j = 1,2,\ldots,n`, by choosing :math:`a` such that a prescribed relative error should be achieved.
The method is modified slightly if :math:`t = 0.0` so that an estimate of :math:`f\left(0.0\right)` can be obtained when it has a finite value. :math:`\tau` is calculated as :math:`t_{\mathrm{fac}}\times \mathrm{max}\left(0.01, t_j\right)`, where :math:`t_{\mathrm{fac}} > 0.5`.
You specify :math:`t_{\mathrm{fac}}` and :math:`\alpha_b`, an upper bound on the exponential order :math:`\alpha` of the inverse function :math:`f\left(t\right)`. :math:`\alpha` has two alternative interpretations:
(i) :math:`\alpha` is the smallest number such that
.. math::
\left\lvert f\left(t\right)\right\rvert \leq m\times \mathrm{exp}\left(\alpha t\right)
for large :math:`t`,
(#) :math:`\alpha` is the real part of the singularity of :math:`F\left(p\right)` with largest real part.
The method depends critically on the value of :math:`\alpha`.
See `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06laf.html#fcomments>`__ for further details.
The function calculates at least two different values of the argument :math:`a`, such that :math:`a > \alpha_b`, in an attempt to achieve the requested relative error and provide error estimates.
The values of :math:`t_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,n`, must be supplied in monotonically increasing order.
The function calculates the values of the inverse function :math:`f\left(t_j\right)` in decreasing order of :math:`j`.
.. _c06la-py2-py-references:
**References**
Crump, K S, 1976, `Numerical inversion of Laplace transforms using a Fourier series approximation`, J. Assoc. Comput. Mach. (23), 89--96
Durbin, F, 1974, `Numerical inversion of Laplace transforms: An efficient improvement to Dubner and Abate's method`, Comput. J. (17), 371--376
Wynn, P, 1956, `On a device for computing the` :math:`e_m\left(S_n\right)` `transformation`, Math. Tables Aids Comput. (10), 91--96
"""
raise NotImplementedError
[docs]def invlaplace_weeks(f, sigma0, sigma, b, epstol, mmax=1024, data=None):
r"""
``invlaplace_weeks`` computes the inverse Laplace transform :math:`f\left(t\right)` of a user-supplied function :math:`F\left(s\right)`, defined for complex :math:`s`.
The function uses a modification of Weeks' method which is suitable when :math:`f\left(t\right)` has continuous derivatives of all orders.
The function returns the coefficients of an expansion which approximates :math:`f\left(t\right)` and can be evaluated for given values of :math:`t` by subsequent calls of :meth:`invlaplace_weeks_eval`.
.. _c06lb-py2-py-doc:
For full information please refer to the NAG Library document for c06lb
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06lbf.html
.. _c06lb-py2-py-parameters:
**Parameters**
**f** : callable retval = f(s, data=None)
:math:`\mathrm{f}` must return the value of the Laplace transform function :math:`F\left(s\right)` for a given complex value of :math:`s`.
**Parameters**
**s** : complex
The value of :math:`s` for which :math:`F\left(s\right)` must be evaluated. The real part of :math:`\mathrm{s}` is greater than :math:`\sigma_0`.
**data** : arbitrary, optional, modifiable in place
User-communication data for callback functions.
**Returns**
**retval** : complex
The value of the Laplace transform function :math:`F\left(s\right)` for the given complex value, :math:`s`.
**sigma0** : float
The abscissa of convergence of the Laplace transform, :math:`\sigma_0`.
**sigma** : float
The parameter :math:`\sigma` of the Laguerre expansion. If on entry :math:`\mathrm{sigma}\leq \sigma_0`, :math:`\mathrm{sigma}` is reset to :math:`\sigma_0+0.7`.
**b** : float
The parameter :math:`b` of the Laguerre expansion. If on entry :math:`\mathrm{b} < 2\left(\sigma -\sigma_0\right)`, :math:`\mathrm{b}` is reset to :math:`2.5\left(\sigma -\sigma_0\right)`.
**epstol** : float
The required relative pseudo-accuracy, that is, an upper bound on :math:`\left\lvert f\left(t\right)-\tilde{f}\left(t\right)\right\rvert e^{{-\sigma t}}`.
**mmax** : int, optional
An upper bound on the number of Laguerre expansion coefficients to be computed. The number of coefficients actually computed is always a power of :math:`2`, so :math:`\mathrm{mmax}` should be a power of :math:`2`; if :math:`\mathrm{mmax}` is not a power of :math:`2` then the maximum number of coefficients calculated will be the largest power of :math:`2` less than :math:`\mathrm{mmax}`.
**data** : arbitrary, optional
User-communication data for callback functions.
**Returns**
**sigma** : float
The value actually used for :math:`\sigma`, as just described.
**b** : float
The value actually used for :math:`b`, as just described.
**m** : int
The number of Laguerre expansion coefficients actually computed. The number of calls to :math:`\mathrm{f}` is :math:`\mathrm{m}/2+2`.
**acoef** : float, ndarray, shape :math:`\left(\mathrm{mmax}\right)`
The first :math:`\mathrm{m}` elements contain the computed Laguerre expansion coefficients, :math:`a_i`.
**errvec** : float, ndarray, shape :math:`\left(8\right)`
An :math:`8`-component vector of diagnostic information.
:math:`\mathrm{errvec}[0]`
Overall estimate of the pseudo-error :math:`\left\lvert f\left(t\right)-\tilde{f}\left(t\right)\right\rvert e^{{-\sigma t}} = \mathrm{errvec}[1]+\mathrm{errvec}[2]+\mathrm{errvec}[3]`.
:math:`\mathrm{errvec}[1]`
Estimate of the discretization pseudo-error.
:math:`\mathrm{errvec}[2]`
Estimate of the truncation pseudo-error.
:math:`\mathrm{errvec}[3]`
Estimate of the condition pseudo-error on the basis of minimal noise levels in function values.
:math:`\mathrm{errvec}[4]`
:math:`K`, coefficient of a heuristic decay function for the expansion coefficients.
:math:`\mathrm{errvec}[5]`
:math:`R`, base of the decay function for the expansion coefficients.
:math:`\mathrm{errvec}[6]`
Logarithm of the largest expansion coefficient.
:math:`\mathrm{errvec}[7]`
Logarithm of the smallest nonzero expansion coefficient.
The values :math:`K` and :math:`R` returned in :math:`\mathrm{errvec}[4]` and :math:`\mathrm{errvec}[5]` define a decay function :math:`KR^{{-i}}` constructed by the function for the purposes of error estimation.
It satisfies
.. math::
\left\lvert a_i\right\rvert < KR^{{-i}}\text{, }i = 1,2,\ldots,m\text{.}
.. _c06lb-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{mmax} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{mmax} \geq 8`.
(`errno` :math:`4`)
The decay rate of the coefficients is too small. Increasing :math:`\mathrm{mmax}` may help. :math:`\mathrm{mmax} = \langle\mathit{\boldsymbol{value}}\rangle`. Pseudo-error estimate :math:`\mathrm{errvec}[0] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{epstol} = \langle\mathit{\boldsymbol{value}}\rangle`.
(`errno` :math:`5`)
The decay rate of the coefficients is too small and round-off error is such that the required accuracy cannot be obtained. Increasing :math:`\mathrm{mmax}` or :math:`\mathrm{epstol}` may help. :math:`\mathrm{mmax} = \langle\mathit{\boldsymbol{value}}\rangle`. Pseudo-error estimate :math:`\mathrm{errvec}[0] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{epstol} = \langle\mathit{\boldsymbol{value}}\rangle`.
**Warns**
**NagAlgorithmicWarning**
(`errno` :math:`2`)
The estimate of the pseudo-error is slightly larger than :math:`\mathrm{epstol}`. Pseudo-error estimate :math:`\mathrm{errvec}[0] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{epstol} = \langle\mathit{\boldsymbol{value}}\rangle`.
(`errno` :math:`3`)
The round-off error level is larger than :math:`\mathrm{epstol}`. Increasing :math:`\mathrm{epstol}` may help. Pseudo-error estimate :math:`\mathrm{errvec}[0] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{epstol} = \langle\mathit{\boldsymbol{value}}\rangle`.
(`errno` :math:`6`)
Error bounds cannot be predicted. Check :math:`\mathrm{sigma0}`. :math:`\mathrm{sigma0} = \langle\mathit{\boldsymbol{value}}\rangle`.
.. _c06lb-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given a function :math:`f\left(t\right)` of a real variable :math:`t`, its Laplace transform :math:`F\left(s\right)` is a function of a complex variable :math:`s`, defined by
.. math::
F\left(s\right) = \int_0^{\infty }e^{{-st}}f\left(t\right){dt}\text{, }\quad \mathrm{Re}\left(s\right) > \sigma_0\text{.}
Then :math:`f\left(t\right)` is the inverse Laplace transform of :math:`F\left(s\right)`.
The value :math:`\sigma_0` is referred to as the abscissa of convergence of the Laplace transform; it is the rightmost real part of the singularities of :math:`F\left(s\right)`.
``invlaplace_weeks``, along with its companion :meth:`invlaplace_weeks_eval`, attempts to solve the following problem:
given a function :math:`F\left(s\right)`, compute values of its inverse Laplace transform :math:`f\left(t\right)` for specified values of :math:`t`.
The method is a modification of Weeks' method (see Garbow `et al.` (1988a)), which approximates :math:`f\left(t\right)` by a truncated Laguerre expansion:
.. math::
\tilde{f}\left(t\right) = e^{{\sigma t}}\sum_{{i = 0}}^{{m-1}}a_ie^{{{-bt}/2}}L_i\left(bt\right)\text{, }\quad \sigma > \sigma_0\text{, }\quad b > 0
where :math:`L_i\left(x\right)` is the Laguerre polynomial of degree :math:`i`.
This function computes the coefficients :math:`a_i` of the above Laguerre expansion; the expansion can then be evaluated for specified :math:`t` by calling :meth:`invlaplace_weeks_eval`.
You must supply the value of :math:`\sigma_0`, and also suitable values for :math:`\sigma` and :math:`b`: see `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06lbf.html#fcomments>`__ for guidance.
The method is only suitable when :math:`f\left(t\right)` has continuous derivatives of all orders.
For such functions the approximation :math:`\tilde{f}\left(t\right)` is usually good and inexpensive.
The function will fail with an error exit if the method is not suitable for the supplied function :math:`F\left(s\right)`.
The function is designed to satisfy an accuracy criterion of the form:
.. math::
\left\lvert \frac{{f\left(t\right)-\tilde{f}\left(t\right)}}{e^{{\sigma t}}}\right\rvert < \epsilon_{\textit{tol}}\text{, for all }t
where :math:`\epsilon_{\textit{tol}}` is a user-supplied bound.
The error measure on the left-hand side is referred to as the **pseudo-relative** **error**, or **pseudo-error** for short.
Note that if :math:`\sigma > 0` and :math:`t` is large, the absolute error in :math:`\tilde{f}\left(t\right)` may be very large.
``invlaplace_weeks`` is derived from the function MODUL1 in Garbow `et al.` (1988a).
.. _c06lb-py2-py-references:
**References**
Garbow, B S, Giunta, G, Lyness, J N and Murli, A, 1988, `Software for an implementation of Weeks' method for the inverse laplace transform problem`, ACM Trans. Math. Software (14), 163--170
Garbow, B S, Giunta, G, Lyness, J N and Murli, A, 1988, `Algorithm 662: A Fortran software package for the numerical inversion of the Laplace transform based on Weeks' method`, ACM Trans. Math. Software (14), 171--176
"""
raise NotImplementedError
[docs]def invlaplace_weeks_eval(t, sigma, b, acoef, errvec):
r"""
``invlaplace_weeks_eval`` evaluates an inverse Laplace transform at a given point, using the expansion coefficients computed by :meth:`invlaplace_weeks`.
.. _c06lc-py2-py-doc:
For full information please refer to the NAG Library document for c06lc
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06lcf.html
.. _c06lc-py2-py-parameters:
**Parameters**
**t** : float
The value :math:`t` for which the inverse Laplace transform :math:`f\left(t\right)` must be evaluated.
**sigma** : float
:math:`\mathrm{sigma}` must be unchanged from the previous call of :meth:`invlaplace_weeks`
**b** : float
:math:`\mathrm{b}` must be unchanged from the previous call of :meth:`invlaplace_weeks`
**acoef** : float, array-like, shape :math:`\left(m\right)`
:math:`\mathrm{acoef}` must be unchanged from the previous call of :meth:`invlaplace_weeks`
**errvec** : float, array-like, shape :math:`\left(8\right)`
:math:`\mathrm{errvec}` must be unchanged from the previous call of :meth:`invlaplace_weeks`
**Returns**
**finv** : float
The approximation to the inverse Laplace transform at :math:`t`.
.. _c06lc-py2-py-errors:
**Warns**
**NagAlgorithmicWarning**
(`errno` :math:`1`)
The approximation to :math:`f\left(t\right)` is too large to be representable.
(`errno` :math:`2`)
The approximation to :math:`f\left(t\right)` is too small to be representable.
.. _c06lc-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``invlaplace_weeks_eval`` is designed to be used following a call to :meth:`invlaplace_weeks`, which computes an inverse Laplace transform by representing it as a Laguerre expansion of the form:
.. math::
\tilde{f}\left(t\right) = e^{{\sigma t}}\sum_{{i = 0}}^{{m-1}}a_ie^{{-{bt}/2}}L_i\left(bt\right)\text{, }\quad \sigma > \sigma_O\text{, }\quad b > 0
where :math:`L_i\left(x\right)` is the Laguerre polynomial of degree :math:`i`.
This function simply evaluates the above expansion for a specified value of :math:`t`.
``invlaplace_weeks_eval`` is derived from the function MODUL2 in Garbow `et al.` (1988)
.. _c06lc-py2-py-references:
**References**
Garbow, B S, Giunta, G, Lyness, J N and Murli, A, 1988, `Algorithm 662: A Fortran software package for the numerical inversion of the Laplace transform based on Weeks' method`, ACM Trans. Math. Software (14), 171--176
"""
raise NotImplementedError
[docs]def fft_realherm_1d(direct, x):
r"""
``fft_realherm_1d`` calculates the discrete Fourier transform of a sequence of :math:`n` real data values or of a Hermitian sequence of :math:`n` complex data values stored in compact form in a `float` array.
.. _c06pa-py2-py-doc:
For full information please refer to the NAG Library document for c06pa
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06paf.html
.. _c06pa-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06pa-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**x** : float, array-like, shape :math:`\left(n+2\right)`
If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n+1\right)` in the function from which ``fft_realherm_1d`` is called:
if :math:`\mathrm{direct} = \texttt{'F'}`, :math:`\mathrm{x}[\textit{j}]` must contain :math:`x_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`;
if :math:`\mathrm{direct} = \texttt{'B'}`, :math:`\mathrm{x}[2\times \textit{k}]` and :math:`\mathrm{x}[2\times {\textit{k}+1}]` must contain the real and imaginary parts respectively of :math:`z_{\textit{k}}`, for :math:`\textit{k} = 0,1,\ldots,{n/2}`. (Note that for the sequence :math:`z_k` to be Hermitian, the imaginary part of :math:`z_0`, and of :math:`z_{{n/2}}` for :math:`n` even, must be zero.)
**Returns**
**x** : float, ndarray, shape :math:`\left(n+2\right)`
If :math:`\mathrm{direct} = \texttt{'F'}` and :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n+1\right)`, :math:`\mathrm{x}[2\times \textit{k}]` and :math:`\mathrm{x}[2\times {\textit{k}+1}]` will contain the real and imaginary parts respectively of :math:`\hat{z}_{\textit{k}}`, for :math:`\textit{k} = 0,1,\ldots,{n/2}`;
if :math:`\mathrm{direct} = \texttt{'B'}` and :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n+1\right)`, :math:`\mathrm{x}[\textit{j}]` will contain :math:`\hat{x}_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
.. _c06pa-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`3`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pa-py2-py-notes:
**Notes**
Given a sequence of :math:`n` real data values :math:`x_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,n-1`, ``fft_realherm_1d`` calculates their discrete Fourier transform (in the **forward** direction) defined by
.. math::
\hat{z}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}x_j\times \mathrm{exp}\left({-i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{.}
The transformed values :math:`\hat{z}_k` are complex, but they form a Hermitian sequence (i.e., :math:`\hat{z}_{{n-k}}` is the complex conjugate of :math:`\hat{z}_k`), so they are completely determined by :math:`n` real numbers (since :math:`\hat{z}_0` is real, as is :math:`\hat{z}_{{n/2}}` for :math:`n` even).
Alternatively, given a Hermitian sequence of :math:`n` complex data values :math:`z_j`, this function calculates their inverse (**backward**) discrete Fourier transform defined by
.. math::
\hat{x}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j\times \mathrm{exp}\left(i\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{.}
The transformed values :math:`\hat{x}_k` are real.
(Note the scale factor of :math:`\frac{1}{\sqrt{n}}` in the above definitions.)
A call of ``fft_realherm_1d`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
``fft_realherm_1d`` uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983).
The same functionality is available using the forward and backward transform function pair: :meth:`fft_real_2d` and :meth:`fft_hermitian_2d` on setting :math:`{\textit{n}} = 1`.
This pair use a different storage solution; real data is stored in a `float` array, while Hermitian data (the first unconjugated half) is stored in a `complex` array.
.. _c06pa-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def fft_complex_1d(direct, x):
r"""
``fft_complex_1d`` calculates the discrete Fourier transform of a sequence of :math:`n` complex data values (using complex data type).
.. _c06pc-py2-py-doc:
For full information please refer to the NAG Library document for c06pc
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pcf.html
.. _c06pc-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06pc-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**x** : complex, array-like, shape :math:`\left(n\right)`
If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_complex_1d`` is called, :math:`\mathrm{x}[\textit{j}]` must contain :math:`z_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
**Returns**
**x** : complex, ndarray, shape :math:`\left(n\right)`
The components of the discrete Fourier transform. If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``fft_complex_1d`` is called, :math:`\hat{z}_k` is contained in :math:`\mathrm{x}[k]`, for :math:`0\leq k\leq {n-1}`.
.. _c06pc-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pc-py2-py-notes:
**Notes**
Given a sequence of :math:`n` complex data values :math:`z_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_complex_1d`` calculates their (**forward** or **backward**) discrete Fourier transform (DFT) defined by
.. math::
\hat{z}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j\times \mathrm{exp}\left({\pm i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{.}
(Note the scale factor of :math:`\frac{1}{\sqrt{n}}` in this definition.) The minus sign is taken in the argument of the exponential within the summation when the forward transform is required, and the plus sign is taken when the backward transform is required.
A call of ``fft_complex_1d`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
``fft_complex_1d`` uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983).
If :math:`n` is a large prime number or if :math:`n` contains large prime factors, then the Fourier transform is performed using Bluestein's algorithm (see Bluestein (1968)), which expresses the DFT as a convolution that in turn can be efficiently computed using FFTs of highly composite sizes.
.. _c06pc-py2-py-references:
**References**
Bluestein, L I, 1968, `A linear filtering approach to the computation of the discrete Fourier transform`, Northeast Electronics Research and Engineering Meeting Record 10, 218--219
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def fft_complex_multid_1d(direct, l, nd, x):
r"""
``fft_complex_multid_1d`` computes the discrete Fourier transform of one variable in a multivariate sequence of complex data values.
.. _c06pf-py2-py-doc:
For full information please refer to the NAG Library document for c06pf
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pff.html
.. _c06pf-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06pf-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**l** : int
:math:`l`, the index of the variable (or dimension) on which the discrete Fourier transform is to be performed.
**nd** : int, array-like, shape :math:`\left(\textit{ndim}\right)`
The elements of :math:`\mathrm{nd}` must contain the dimensions of the :math:`\textit{ndim}` variables; that is, :math:`\mathrm{nd}[i-1]` must contain the dimension of the :math:`i`\ th variable.
**x** : complex, array-like, shape :math:`\left(n\right)`
The complex data values. Data values are stored in :math:`\mathrm{x}` using column-major ordering for storing multidimensional arrays; that is, :math:`z_{{j_1j_2 \cdots j_m}}` is stored in :math:`\mathrm{x}[ 1 + j_1 + n_1 j_2 + n_1 n_2 j_3 + \cdots ]`.
**Returns**
**x** : complex, ndarray, shape :math:`\left(n\right)`
The corresponding elements of the computed transform.
.. _c06pf-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\textit{ndim} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\textit{ndim} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{l} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`1\leq \mathrm{l}\leq \textit{ndim}`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{nd}[\langle\mathit{\boldsymbol{value}}\rangle] = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{nd}[i-1] \geq 1`, for all :math:`i`.
(`errno` :math:`5`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`, product of :math:`\mathrm{nd}` elements is :math:`\langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\textit{n}` must equal the product of the dimensions held in array :math:`\mathrm{nd}`.
(`errno` :math:`8`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pf-py2-py-notes:
**Notes**
``fft_complex_multid_1d`` computes the discrete Fourier transform of one variable (the :math:`l`\ th say) in a multivariate sequence of complex data values :math:`z_{{j_1j_2 \cdots j_m}}`, where :math:`j_1 = 0,1,\ldots,{n_1-1}\text{, }\quad j_2 = 0,1,\ldots,{n_2-1}`, and so on.
Thus the individual dimensions are :math:`n_1,n_2,\ldots,n_m`, and the total number of data values is :math:`n = n_1\times n_2\times \cdots \times n_m`.
The function computes :math:`n/n_l` one-dimensional transforms defined by
.. math::
\hat{z}_{{j_1\ldots k_l\ldots j_m}} = \frac{1}{\sqrt{n_l}}\sum_{{j_l = 0}}^{{n_l-1}}z_{{j_1\ldots j_l\ldots j_m}}\times \mathrm{exp}\left(\pm \frac{{2\pi ij_lk_l}}{n_l}\right)\text{,}
where :math:`k_l = 0,1,\ldots,n_l-1`.
The plus or minus sign in the argument of the exponential terms in the above definition determine the direction of the transform: a minus sign defines the **forward** direction and a plus sign defines the **backward** direction.
(Note the scale factor of :math:`\frac{1}{\sqrt{n_l}}` in this definition.)
A call of ``fft_complex_multid_1d`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
The data values must be supplied in a one-dimensional complex array using column-major storage ordering of multidimensional data (i.e., with the first subscript :math:`j_1` varying most rapidly).
This function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983).
.. _c06pf-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def fft_complex_multid(direct, nd, x):
r"""
``fft_complex_multid`` computes the multidimensional discrete Fourier transform of a multivariate sequence of complex data values.
.. _c06pj-py2-py-doc:
For full information please refer to the NAG Library document for c06pj
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pjf.html
.. _c06pj-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06pj-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**nd** : int, array-like, shape :math:`\left(\textit{ndim}\right)`
The elements of :math:`\mathrm{nd}` must contain the dimensions of the :math:`\textit{ndim}` variables; that is, :math:`\mathrm{nd}[i-1]` must contain the dimension of the :math:`i`\ th variable.
**x** : complex, array-like, shape :math:`\left(n\right)`
The complex data values. Data values are stored in :math:`\mathrm{x}` using column-major ordering for storing multidimensional arrays; that is, :math:`z_{{j_1j_2 \cdots j_m}}` is stored in :math:`\mathrm{x}[ 1 + j_1 + n_1 j_2 + n_1 n_2 j_3 + \cdots ]`.
**Returns**
**x** : complex, ndarray, shape :math:`\left(n\right)`
The corresponding elements of the computed transform.
.. _c06pj-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\textit{ndim} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\textit{ndim} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{nd}[\langle\mathit{\boldsymbol{value}}\rangle] = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{nd}[i-1] \geq 1`, for all :math:`i`.
(`errno` :math:`4`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`, product of :math:`\mathrm{nd}` elements is :math:`\langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\textit{n}` must equal the product of the dimensions held in array :math:`\mathrm{nd}`.
(`errno` :math:`7`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pj-py2-py-notes:
**Notes**
``fft_complex_multid`` computes the multidimensional discrete Fourier transform of a multidimensional sequence of complex data values :math:`z_{{j_1j_2\ldots j_m}}`, where :math:`j_1 = 0,1,\ldots,n_1-1\text{, }\quad j_2 = 0,1,\ldots,n_2-1`, and so on.
Thus the individual dimensions are :math:`n_1,n_2,\ldots,n_m`, and the total number of data values is :math:`n = n_1\times n_2\times \cdots \times n_m`.
The discrete Fourier transform is here defined (e.g., for :math:`m = 2`) by:
.. math::
\hat{z}_{{k_1,k_2}} = \frac{1}{\sqrt{n}}\sum_{{j_1 = 0}}^{{n_1-1}}\sum_{{j_2 = 0}}^{{n_2-1}}z_{{j_1j_2}}\times \mathrm{exp}\left({\pm 2\pi i}\left(\frac{{j_1k_1}}{n_1}+\frac{{j_2k_2}}{n_2}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,n_1-1` and :math:`k_2 = 0,1,\ldots,n_2-1`.
The plus or minus sign in the argument of the exponential terms in the above definition determine the direction of the transform: a minus sign defines the **forward** direction and a plus sign defines the **backward** direction.
The extension to higher dimensions is obvious. (Note the scale factor of :math:`\frac{1}{\sqrt{n}}` in this definition.)
A call of ``fft_complex_multid`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
The data values must be supplied in a one-dimensional array using column-major storage ordering of multidimensional data (i.e., with the first subscript :math:`j_1` varying most rapidly).
This function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983).
.. _c06pj-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def convcorr_complex(job, x, y):
r"""
``convcorr_complex`` calculates the circular convolution or correlation of two complex vectors of period :math:`n`.
.. _c06pk-py2-py-doc:
For full information please refer to the NAG Library document for c06pk
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pkf.html
.. _c06pk-py2-py-parameters:
**Parameters**
**job** : int
The computation to be performed:
:math:`\mathrm{job} = 1`
:math:`z_k = {\sum_{{j = 0}}^{{n-1}}x_jy_{{k-j}}}` (convolution);
:math:`\mathrm{job} = 2`
:math:`w_k = {\sum_{{j = 0}}^{{n-1}}\bar{x}_jy_{{k+j}}}` (correlation).
**x** : complex, array-like, shape :math:`\left(n\right)`
The elements of one period of the vector :math:`x`. If :math:`\mathrm{x}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``convcorr_complex`` is called, :math:`\mathrm{x}[\textit{j}]` must contain :math:`x_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
**y** : complex, array-like, shape :math:`\left(n\right)`
The elements of one period of the vector :math:`y`. If :math:`\mathrm{y}` is declared with bounds :math:`\left(0:n-1\right)` in the function from which ``convcorr_complex`` is called, :math:`\mathrm{y}[\textit{j}-1]` must contain :math:`y_{\textit{j}}`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`.
**Returns**
**x** : complex, ndarray, shape :math:`\left(n\right)`
The corresponding elements of the discrete convolution or correlation.
**y** : complex, ndarray, shape :math:`\left(n\right)`
The discrete Fourier transform of the convolution or correlation returned in the array :math:`\mathrm{x}`.
.. _c06pk-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`2`)
On entry, argument :math:`\mathrm{job}` had an illegal value.
.. _c06pk-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``convcorr_complex`` computes:
if :math:`\mathrm{job} = 1`, the discrete **convolution** of :math:`x` and :math:`y`, defined by
.. math::
z_k = \sum_{{j = 0}}^{{n-1}}x_jy_{{k-j}} = \sum_{{j = 0}}^{{n-1}}x_{{k-j}}y_j\text{;}
if :math:`\mathrm{job} = 2`, the discrete **correlation** of :math:`x` and :math:`y` defined by
.. math::
w_k = \sum_{{j = 0}}^{{n-1}}\bar{x}_jy_{{k+j}}\text{.}
Here :math:`x` and :math:`y` are complex vectors, assumed to be periodic, with period :math:`n`, i.e., :math:`x_j = x_{{j\pm n}} = x_{{j\pm 2n}} = \cdots \text{}`; :math:`z` and :math:`w` are then also periodic with period :math:`n`.
**Note:** this usage of the terms 'convolution' and 'correlation' is taken from Brigham (1974).
The term 'convolution' is sometimes used to denote both these computations.
If :math:`\hat{x}`, :math:`\hat{y}`, :math:`\hat{z}` and :math:`\hat{w}` are the discrete Fourier transforms of these sequences, and :math:`\tilde{x}` is the inverse discrete Fourier transform of the sequence :math:`x_j`, i.e.,
.. math::
\hat{x}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}x_j\times \mathrm{exp}\left({-i}\frac{{2\pi jk}}{n}\right)\text{, etc.,}
and
.. math::
\tilde{x}_k = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}x_j\times \mathrm{exp}\left(i\frac{{2\pi jk}}{n}\right)\text{,}
then :math:`\hat{z}_k = \sqrt{n}.\hat{x}_k\hat{y}_k` and :math:`\hat{w}_k = \sqrt{n}.\bar{\hat{x}}_k\hat{y}_k` (the bar denoting complex conjugate).
.. _c06pk-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
"""
raise NotImplementedError
[docs]def fft_realherm_1d_multi_row(direct, m, n, x):
r"""
``fft_realherm_1d_multi_row`` computes the discrete Fourier transforms of :math:`m` sequences, each containing :math:`n` real data values or a Hermitian complex sequence stored in a complex storage format.
.. _c06pp-py2-py-doc:
For full information please refer to the NAG Library document for c06pp
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06ppf.html
.. _c06pp-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06pp-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**m** : int
:math:`m`, the number of sequences to be transformed.
**n** : int
:math:`n`, the number of real or complex values in each sequence.
**x** : float, array-like, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+2\right)\right)`
The data must be stored such that consecutive elements of the same sequence are stored with a stride of :math:`\mathrm{m}` and corresponding elements of different sequences are stored consecutively. An additional two spaces are reserved for each sequence to allow for the pairwise storage of real and imaginary parts in the transformed domain. In other words, if the data values of the :math:`p`\ th sequence to be transformed are denoted by :math:`x_{\textit{j}}^p`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`:
if :math:`\mathrm{direct} = \texttt{'F'}`, :math:`\mathrm{x}[\textit{j}\times \mathrm{m}+\textit{p}-1]` must contain :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`;
if :math:`\mathrm{direct} = \texttt{'B'}`, :math:`\mathrm{x}[2\times \textit{k}\times \mathrm{m}+\textit{p}-1]` and :math:`\mathrm{x}[ \left(2\times {\textit{k}+1}\right) \times {\mathrm{m}+\textit{p}}-1]` must contain the real and imaginary parts respectively of :math:`\hat{z}_k^p`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 0,1,\ldots,{n/2}`. (Note that for the sequence :math:`\hat{z}_k^p` to be Hermitian, the imaginary part of :math:`\hat{z}_0^p`, and of :math:`\hat{z}_{{n/2}}^p` for :math:`n` even, must be zero.)
**Returns**
**x** : float, ndarray, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+2\right)\right)`
If :math:`\mathrm{direct} = \texttt{'F'}` then :math:`\mathrm{x}[2\times \textit{k}\times \mathrm{m}+\textit{p}-1]` and :math:`\mathrm{x}[\left(2\times \textit{k}+1\right)\times \mathrm{m}+\textit{p}-1]` will contain the real and imaginary parts respectively of :math:`\hat{z}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 0,1,\ldots,{n/2}`;
if :math:`\mathrm{direct} = \texttt{'B'}` then :math:`\mathrm{x}[\textit{j}\times \mathrm{m}+\textit{p}-1]` will contain :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`;
.. _c06pp-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pp-py2-py-notes:
**Notes**
Given :math:`m` sequences of :math:`n` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_realherm_1d_multi_row`` simultaneously calculates the Fourier transforms of all the sequences defined by
.. math::
\hat{z}_k^p = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}x_j^p\times \mathrm{exp}\left({-i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{ and }p = 1,2,\ldots,m\text{.}
The transformed values :math:`\hat{z}_k^p` are complex, but for each value of :math:`p` the :math:`\hat{z}_k^p` form a Hermitian sequence (i.e., :math:`\hat{z}_{{n-k}}^p` is the complex conjugate of :math:`\hat{z}_k^p`), so they are completely determined by :math:`mn` real numbers (since :math:`\hat{z}_0^p` is real, as is :math:`\hat{z}_{{n/2}}^p` for :math:`n` even).
Alternatively, given :math:`m` Hermitian sequences of :math:`n` complex data values :math:`z_j^p`, this function simultaneously calculates their inverse (**backward**) discrete Fourier transforms defined by
.. math::
\hat{x}_k^p = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j^p\times \mathrm{exp}\left(i\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{ and }p = 1,2,\ldots,m\text{.}
The transformed values :math:`\hat{x}_k^p` are real.
(Note the scale factor :math:`\frac{1}{\sqrt{n}}` in the above definition.)
A call of ``fft_realherm_1d_multi_row`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06pp-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_realherm_1d_multi_col(direct, n, m, x):
r"""
``fft_realherm_1d_multi_col`` computes the discrete Fourier transforms of :math:`m` sequences, each containing :math:`n` real data values or a Hermitian complex sequence stored column-wise in a complex storage format.
.. _c06pq-py2-py-doc:
For full information please refer to the NAG Library document for c06pq
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pqf.html
.. _c06pq-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06pq-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**n** : int
:math:`n`, the number of real or complex values in each sequence.
**m** : int
:math:`m`, the number of sequences to be transformed.
**x** : float, array-like, shape :math:`\left(\left(\mathrm{n}+2\right)\times \mathrm{m}\right)`
The :math:`m` real or Hermitian data sequences to be transformed.
if :math:`\mathrm{direct} = \texttt{'F'}`, the :math:`m` real data sequences, :math:`x^{\textit{p}} = \left(x_1^{\textit{p}},x_2^{\textit{p}},\ldots,x_n^{\textit{p}}\right)`, for :math:`\textit{p} = 1,2,\ldots,m`, should be stored sequentially in :math:`\mathrm{x}`, with a stride of :math:`n+2` between sequences.
if :math:`\mathrm{direct} = \texttt{'B'}`, the :math:`m` Hermitian data sequences, :math:`\hat{z}^{\textit{p}} = \left(\hat{z}_1^{\textit{p}},\hat{z}_2^{\textit{p}},\ldots,\hat{z}_{{n/2+1}}^{\textit{p}}\right) = \left(\mathrm{Re}\left(\hat{z}_1^{\textit{p}}\right),\mathrm{Im}\left(\hat{z}_1^{\textit{p}}\right),\mathrm{Re}\left(\hat{z}_2^{\textit{p}}\right),\mathrm{Im}\left(\hat{z}_2^{\textit{p}}\right),\ldots,\mathrm{Re}\left(\hat{z}_{{n/2+1}}^{\textit{p}}\right),\mathrm{Im}\left(\hat{z}_{{n/2+1}}^{\textit{p}}\right)\right)`, for :math:`\textit{p} = 1,2,\ldots,m`, should be stored sequentially in :math:`\mathrm{x}`, with a stride of :math:`n+2` between sequences.
In other words:
if :math:`\mathrm{direct} = \texttt{'F'}`, :math:`\mathrm{x}[\left(\textit{p}-1\right)\times \left(\mathrm{n}+2\right)+\textit{j}-1]` must contain :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,n`;
if :math:`\mathrm{direct} = \texttt{'B'}`, :math:`\mathrm{x}[\left(\textit{p}-1\right) \times \left(\mathrm{n}+2\right) + \left(2\times {\textit{k}-1}\right)-1]` and :math:`\mathrm{x}[\left(\textit{p}-1\right) \times \left(\mathrm{n}+2\right) + \left(2\times \textit{k}\right)-1]` must contain the real and imaginary parts respectively of :math:`\hat{z}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 1,2,\ldots,{n/2+1}`. (Note that for the sequence :math:`\hat{z}_k^p` to be Hermitian, the imaginary part of :math:`\hat{z}_1^p`, and of :math:`\hat{z}_{{n/2+1}}^p` for :math:`n` even, must be zero.)
**Returns**
**x** : float, ndarray, shape :math:`\left(\left(\mathrm{n}+2\right)\times \mathrm{m}\right)`
If :math:`\mathrm{direct} = \texttt{'F'}` then the :math:`m` sequences, :math:`\hat{z}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m` stored as described on entry for :math:`\mathrm{direct} = \texttt{'B'}`
if :math:`\mathrm{direct} = \texttt{'B'}` then the :math:`m` sequences, :math:`x^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m` stored as described on entry for :math:`\mathrm{direct} = \texttt{'F'}`
.. _c06pq-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pq-py2-py-notes:
**Notes**
Given :math:`m` sequences of :math:`n` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_realherm_1d_multi_col`` simultaneously calculates the Fourier transforms of all the sequences defined by
.. math::
\hat{z}_k^{\textit{p}} = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}x_j^p\times \mathrm{exp}\left({-i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{ and }p = 1,2,\ldots,m\text{.}
The transformed values :math:`\hat{z}_k^p` are complex, but for each value of :math:`p` the :math:`\hat{z}_k^p` form a Hermitian sequence (i.e., :math:`\hat{z}_{{n-k}}^p` is the complex conjugate of :math:`\hat{z}_k^p`), so they are completely determined by :math:`mn` real numbers (since :math:`\hat{z}_0^p` is real, as is :math:`\hat{z}_{{n/2}}^p` for :math:`n` even).
Alternatively, given :math:`m` Hermitian sequences of :math:`n` complex data values :math:`z_j^p`, this function simultaneously calculates their inverse (**backward**) discrete Fourier transforms defined by
.. math::
\hat{x}_k^p = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j^p\times \mathrm{exp}\left(i\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{ and }p = 1,2,\ldots,m\text{.}
The transformed values :math:`\hat{x}_k^p` are real.
(Note the scale factor :math:`\frac{1}{\sqrt{n}}` in the above definition.)
A call of ``fft_realherm_1d_multi_col`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06pq-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_complex_1d_multi_row(direct, m, n, x):
r"""
``fft_complex_1d_multi_row`` computes the discrete Fourier transforms of :math:`m` sequences, each containing :math:`n` complex data values.
.. _c06pr-py2-py-doc:
For full information please refer to the NAG Library document for c06pr
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06prf.html
.. _c06pr-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06pr-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**m** : int
:math:`m`, the number of sequences to be transformed.
**n** : int
:math:`n`, the number of complex values in each sequence.
**x** : complex, array-like, shape :math:`\left(\mathrm{m}\times \mathrm{n}\right)`
The complex data must be stored in :math:`\mathrm{x}` as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {0:\mathrm{n}-1}\right)`; each of the :math:`m` sequences is stored in a **row** of each array. In other words, if the elements of the :math:`p`\ th sequence to be transformed are denoted by :math:`z_{\textit{j}}^p`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, :math:`\mathrm{x}[j\times \mathrm{m}+p-1]` must contain :math:`z_j^p`.
**Returns**
**x** : complex, ndarray, shape :math:`\left(\mathrm{m}\times \mathrm{n}\right)`
Is overwritten by the complex transforms.
.. _c06pr-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`5`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pr-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given :math:`m` sequences of :math:`n` complex data values :math:`z_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_complex_1d_multi_row`` simultaneously calculates the (**forward** or **backward**) discrete Fourier transforms of all the sequences defined by
.. math::
\hat{z}_k^p = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j^p\times \mathrm{exp}\left({\pm i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{ and }p = 1,2,\ldots,m\text{.}
(Note the scale factor :math:`\frac{1}{\sqrt{n}}` in this definition.) The minus sign is taken in the argument of the exponential within the summation when the forward transform is required, and the plus sign is taken when the backward transform is required.
A call of ``fft_complex_1d_multi_row`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983).
Special code is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06pr-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def fft_complex_1d_multi_col(direct, n, m, x):
r"""
``fft_complex_1d_multi_col`` computes the discrete Fourier transforms of :math:`m` sequences, stored as columns of an array, each containing :math:`n` complex data values.
.. _c06ps-py2-py-doc:
For full information please refer to the NAG Library document for c06ps
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06psf.html
.. _c06ps-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06ps-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**n** : int
:math:`n`, the number of complex values in each sequence.
**m** : int
:math:`m`, the number of sequences to be transformed.
**x** : complex, array-like, shape :math:`\left(\mathrm{n}\times \mathrm{m}\right)`
The complex data values :math:`z_{\textit{j}}^p` stored in :math:`\mathrm{x}[\left(\textit{p}-1\right)\times \mathrm{n}+\textit{j}]`, for :math:`\textit{p} = 1,2,\ldots,\mathrm{m}`, for :math:`\textit{j} = 0,1,\ldots,\mathrm{n}-1`.
**Returns**
**x** : complex, ndarray, shape :math:`\left(\mathrm{n}\times \mathrm{m}\right)`
Is overwritten by the complex transforms.
.. _c06ps-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`5`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06ps-py2-py-notes:
**Notes**
Given :math:`m` sequences of :math:`n` complex data values :math:`z_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_complex_1d_multi_col`` simultaneously calculates the (**forward** or **backward**) discrete Fourier transforms of all the sequences defined by
.. math::
\hat{z}_k^p = \frac{1}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}z_j^p\times \mathrm{exp}\left({\pm i}\frac{{2\pi jk}}{n}\right)\text{, }\quad k = 0,1,\ldots,{n-1}\text{ and }p = 1,2,\ldots,m\text{.}
(Note the scale factor :math:`\frac{1}{\sqrt{n}}` in this definition.) The minus sign is taken in the argument of the exponential within the summation when the forward transform is required, and the plus sign is taken when the backward transform is required.
A call of ``fft_complex_1d_multi_col`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983).
Special code is provided for the factors :math:`2`, :math:`3` and :math:`5`.
.. _c06ps-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def fft_complex_2d(direct, m, n, x):
r"""
``fft_complex_2d`` computes the two-dimensional discrete Fourier transform of a bivariate sequence of complex data values (using complex data type).
.. _c06pu-py2-py-doc:
For full information please refer to the NAG Library document for c06pu
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06puf.html
.. _c06pu-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06pu-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**m** : int
:math:`m`, the first dimension of the transform.
**n** : int
:math:`n`, the second dimension of the transform.
**x** : complex, array-like, shape :math:`\left(\mathrm{m}\times \mathrm{n}\right)`
The complex data values. :math:`\mathrm{x}[\mathrm{m}\times \textit{j}_2+\textit{j}_1]` must contain :math:`z_{{\textit{j}_1\textit{j}_2}}`, for :math:`\textit{j}_2 = 0,\ldots,\mathrm{n}-1`, for :math:`\textit{j}_1 = 0,\ldots,\mathrm{m}-1`.
**Returns**
**x** : complex, ndarray, shape :math:`\left(\mathrm{m}\times \mathrm{n}\right)`
The corresponding elements of the computed transform.
.. _c06pu-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`6`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pu-py2-py-notes:
**Notes**
``fft_complex_2d`` computes the two-dimensional discrete Fourier transform of a bivariate sequence of complex data values :math:`z_{{\textit{j}_1\textit{j}_2}}`, for :math:`\textit{j}_2 = 0,1,\ldots,{n-1}`, for :math:`\textit{j}_1 = 0,1,\ldots,{m-1}`.
The discrete Fourier transform is here defined by
.. math::
\hat{z}_{{k_1k_2}} = \frac{1}{\sqrt{mn}}\sum_{{j_1 = 0}}^{{m-1}}\sum_{{j_2 = 0}}^{{n-1}}z_{{j_1j_2}}\times \mathrm{exp}\left({\pm 2\pi i}\left(\frac{{j_1k_1}}{m}+\frac{{j_2k_2}}{n}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,{m-1}` and :math:`k_2 = 0,1,\ldots,{n-1}`.
(Note the scale factor of :math:`\frac{1}{\sqrt{mn}}` in this definition.) The minus sign is taken in the argument of the exponential within the summation when the forward transform is required, and the plus sign is taken when the backward transform is required.
A call of ``fft_complex_2d`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
This function calls :meth:`fft_complex_1d_multi_row` to perform multiple one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm in Brigham (1974).
.. _c06pu-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def fft_real_2d(m, n, x):
r"""
``fft_real_2d`` computes the two-dimensional discrete Fourier transform of a bivariate sequence of real data values.
.. _c06pv-py2-py-doc:
For full information please refer to the NAG Library document for c06pv
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pvf.html
.. _c06pv-py2-py-parameters:
**Parameters**
**m** : int
:math:`m`, the first dimension of the transform.
**n** : int
:math:`n`, the second dimension of the transform.
**x** : float, array-like, shape :math:`\left(\mathrm{m}\times \mathrm{n}\right)`
The real input dataset :math:`x`, where :math:`x_{{\textit{j}_1\textit{j}_2}}` is stored in :math:`\mathrm{x}[ \textit{j}_2 \times m + \textit{j}_1]`, for :math:`\textit{j}_2 = 0,1,\ldots,n-1`, for :math:`\textit{j}_1 = 0,1,\ldots,m-1`.
**Returns**
**y** : complex, ndarray, shape :math:`\left(\left(\mathrm{m}/2+1\right)\times \mathrm{n}\right)`
The complex output dataset :math:`\hat{z}`, where :math:`\hat{z}_{{\textit{k}_1\textit{k}_2}}` is stored in :math:`\mathrm{y}[ \textit{k}_2 \times \left(m/2+1\right) + k_1]`, for :math:`\textit{k}_2 = 0,1,\ldots,n-1`, for :math:`\textit{k}_1 = 0,1,\ldots,m/2`. Note the first dimension is cut roughly by half to remove the redundant information due to conjugate symmetry.
.. _c06pv-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pv-py2-py-notes:
**Notes**
``fft_real_2d`` computes the two-dimensional discrete Fourier transform of a bivariate sequence of real data values :math:`x_{{\textit{j}_1\textit{j}_2}}`, for :math:`\textit{j}_2 = 0,1,\ldots,n-1`, for :math:`\textit{j}_1 = 0,1,\ldots,m-1`.
The discrete Fourier transform is here defined by
.. math::
\hat{z}_{{k_1k_2}} = \frac{1}{\sqrt{mn}}\sum_{0}^{m-1}{\sum_{0}^{n-1}{x_{{j_1j_2}}}}\times \mathrm{exp}\left({-2\pi i}\left(\frac{{j_1k_1}}{m}+\frac{{j_2k_2}}{n}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,{m-1}` and :math:`k_2 = 0,1,\ldots,{n-1}`. (Note the scale factor of :math:`\frac{1}{\sqrt{mn}}` in this definition.)
The transformed values :math:`\hat{z}_{{k_1k_2}}` are complex.
Because of conjugate symmetry (i.e., :math:`\hat{z}_{{k_1k_2}}` is the complex conjugate of :math:`\hat{z}_{{\left(m-k_1\right)\left(n-k_2\right)}}`), only slightly more than half of the Fourier coefficients need to be stored in the output.
A call of ``fft_real_2d`` followed by a call of :meth:`fft_hermitian_2d` will restore the original data.
This function calls :meth:`fft_realherm_1d_multi_col` and :meth:`fft_complex_1d_multi_row` to perform multiple one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm in Brigham (1974) and Temperton (1983).
.. _c06pv-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_hermitian_2d(m, n, y):
r"""
``fft_hermitian_2d`` computes the two-dimensional inverse discrete Fourier transform of a bivariate Hermitian sequence of complex data values.
.. _c06pw-py2-py-doc:
For full information please refer to the NAG Library document for c06pw
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pwf.html
.. _c06pw-py2-py-parameters:
**Parameters**
**m** : int
:math:`m`, the first dimension of the transform.
**n** : int
:math:`n`, the second dimension of the transform.
**y** : complex, array-like, shape :math:`\left(\left(\mathrm{m}/2+1\right)\times \mathrm{n}\right)`
The Hermitian sequence of complex input dataset :math:`z`, where :math:`z_{{\textit{j}_1\textit{j}_2}}` is stored in :math:`\mathrm{y}[ \textit{j}_2 \times \left(m/2+1\right) + j_1]`, for :math:`\textit{j}_2 = 0,1,\ldots,n-1`, for :math:`\textit{j}_1 = 0,1,\ldots,m/2`.
**Returns**
**x** : float, ndarray, shape :math:`\left(\mathrm{m}\times \mathrm{n}\right)`
The real output dataset :math:`\hat{x}`, where :math:`\hat{x}_{{\textit{k}_1\textit{k}_2}}` is stored in :math:`\mathrm{x}[ \textit{k}_2 \times m + \textit{k}_1]`, for :math:`\textit{k}_2 = 0,1,\ldots,n-1`, for :math:`\textit{k}_1 = 0,1,\ldots,m-1`.
.. _c06pw-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pw-py2-py-notes:
**Notes**
``fft_hermitian_2d`` computes the two-dimensional inverse discrete Fourier transform of a bivariate Hermitian sequence of complex data values :math:`z_{{\textit{j}_1\textit{j}_2}}`, for :math:`\textit{j}_2 = 0,1,\ldots,n-1`, for :math:`\textit{j}_1 = 0,1,\ldots,m-1`.
The discrete Fourier transform is here defined by
.. math::
\hat{x}_{{k_1k_2}} = \frac{1}{\sqrt{mn}}\sum_{0}^{m-1}{\sum_{0}^{n-1}{z_{{j_1j_2}}}}\times \mathrm{exp}\left({2\pi i}\left(\frac{{j_1k_1}}{m}+\frac{{j_2k_2}}{n}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,{m-1}` and :math:`k_2 = 0,1,\ldots,{n-1}`. (Note the scale factor of :math:`\frac{1}{\sqrt{mn}}` in this definition.)
Because the input data satisfies conjugate symmetry (i.e., :math:`z_{{j_1j_2}}` is the complex conjugate of :math:`z_{{\left(m-j_1\right)\left(n-j_2\right)}}`, the transformed values :math:`\hat{x}_{{k_1k_2}}` are real.
A call of :meth:`fft_real_2d` followed by a call of ``fft_hermitian_2d`` will restore the original data.
This function calls :meth:`fft_realherm_1d_multi_col` and :meth:`fft_complex_1d_multi_row` to perform multiple one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm in Brigham (1974) and Temperton (1983).
.. _c06pw-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_complex_3d(direct, n1, n2, n3, x):
r"""
``fft_complex_3d`` computes the three-dimensional discrete Fourier transform of a trivariate sequence of complex data values (using complex data type).
.. _c06px-py2-py-doc:
For full information please refer to the NAG Library document for c06px
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pxf.html
.. _c06px-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06px-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**n1** : int
:math:`n_1`, the first dimension of the transform.
**n2** : int
:math:`n_2`, the second dimension of the transform.
**n3** : int
:math:`n_3`, the third dimension of the transform.
**x** : complex, array-like, shape :math:`\left(\mathrm{n1}\times \mathrm{n2}\times \mathrm{n3}\right)`
The complex data values. Data values are stored in :math:`\mathrm{x}` using column-major ordering for storing multidimensional arrays; that is, :math:`z_{{j_1j_2j_3}}` is stored in :math:`\mathrm{x}[ 1 + j_1 + n_1 j_2 + n_1 n_2 j_3 ]`.
**Returns**
**x** : complex, ndarray, shape :math:`\left(\mathrm{n1}\times \mathrm{n2}\times \mathrm{n3}\right)`
The corresponding elements of the computed transform.
.. _c06px-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{n1} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n1} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n2} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n2} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{n3} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n3} \geq 1`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`8`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06px-py2-py-notes:
**Notes**
``fft_complex_3d`` computes the three-dimensional discrete Fourier transform of a trivariate sequence of complex data values :math:`z_{{j_1j_2j_3}}`, for :math:`\textit{j}_3 = 0,1,\ldots,{n_3-1}`, for :math:`\textit{j}_2 = 0,1,\ldots,{n_2-1}`, for :math:`\textit{j}_1 = 0,1,\ldots,{n_1-1}`.
The discrete Fourier transform is here defined by
.. math::
\hat{z}_{{k_1k_2k_3}} = \frac{1}{\sqrt{n_1n_2n_3}}\sum_{{j_1 = 0}}^{{n_1-1}}\sum_{{j_2 = 0}}^{{n_2-1}}\sum_{{j_3 = 0}}^{{n_3-1}}z_{{j_1j_2j_3}}\times \mathrm{exp}\left({\pm 2\pi i}\left(\frac{{j_1k_1}}{n_1}+\frac{{j_2k_2}}{n_2}+\frac{{j_3k_3}}{n_3}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,{n_1-1}`, :math:`k_2 = 0,1,\ldots,{n_2-1}` and :math:`k_3 = 0,1,\ldots,{n_3-1}`.
(Note the scale factor of :math:`\frac{1}{\sqrt{n_1n_2n_3}}` in this definition.) The minus sign is taken in the argument of the exponential within the summation when the forward transform is required, and the plus sign is taken when the backward transform is required.
A call of ``fft_complex_3d`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
This function performs multiple one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm (see Brigham (1974)).
.. _c06px-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Self-sorting mixed-radix fast Fourier transforms`, J. Comput. Phys. (52), 1--23
"""
raise NotImplementedError
[docs]def fft_real_3d(n1, n2, n3, x):
r"""
``fft_real_3d`` computes the three-dimensional discrete Fourier transform of a trivariate sequence of real data values.
.. _c06py-py2-py-doc:
For full information please refer to the NAG Library document for c06py
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pyf.html
.. _c06py-py2-py-parameters:
**Parameters**
**n1** : int
:math:`n_1`, the first dimension of the transform.
**n2** : int
:math:`n_2`, the second dimension of the transform.
**n3** : int
:math:`n_3`, the third dimension of the transform.
**x** : float, array-like, shape :math:`\left(\mathrm{n1}\times \mathrm{n2}\times \mathrm{n3}\right)`
The real input dataset :math:`x`, where :math:`x_{{\textit{j}_1\textit{j}_2\textit{j}_3}}` is stored in :math:`\mathrm{x}[ \textit{j}_3 \times n_1 n_2 + \textit{j}_2 \times n_1 + \textit{j}_1 +1 -1]`, for :math:`\textit{j}_3 = 0,1,\ldots,n_3-1`, for :math:`\textit{j}_2 = 0,1,\ldots,n_2-1`, for :math:`\textit{j}_1 = 0,1,\ldots,n_1-1`.
**Returns**
**y** : complex, ndarray, shape :math:`\left(\left(\mathrm{n1}/2+1\right)\times \mathrm{n2}\times \mathrm{n3}\right)`
The complex output dataset :math:`\hat{z}`, where :math:`\hat{z}_{{\textit{k}_1\textit{k}_2\textit{k}_3}}` is stored in :math:`\mathrm{y}[ \textit{k}_3 \times \left(n_1/2+1\right) n_2 + \textit{k}_2 \times \left(n_1/2+1\right) + \textit{k}_1 +1 -1]`, for :math:`\textit{k}_3 = 0,1,\ldots,n_3-1`, for :math:`\textit{k}_2 = 0,1,\ldots,n_2-1`, for :math:`\textit{k}_1 = 0,1,\ldots,n_1/2`. Note the first dimension is cut roughly by half to remove the redundant information due to conjugate symmetry.
.. _c06py-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{n1} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n1} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n2} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n2} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{n3} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n3} \geq 1`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06py-py2-py-notes:
**Notes**
``fft_real_3d`` computes the three-dimensional discrete Fourier transform of a trivariate sequence of real data values :math:`x_{{j_1j_2j_3}}`, for :math:`\textit{j}_3 = 0,1,\ldots,n_3-1`, for :math:`\textit{j}_2 = 0,1,\ldots,n_2-1`, for :math:`\textit{j}_1 = 0,1,\ldots,n_1-1`.
The discrete Fourier transform is here defined by
.. math::
\hat{z}_{{k_1k_2k_3}} = \frac{1}{{\sqrt{n_1n_2n_3}}}\sum_{0}^{n_1-1}{\sum_{0}^{n_2-1}{\sum_{0}^{n_3-1}{x_{{j_1j_2j_3}}}}}\times \mathrm{exp}\left(-2πi\left(\frac{{j_1k_1}}{{n_1}}+\frac{{j_2k_2}}{{n_2}}+\frac{{j_3k_3}}{{n_3}}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,n_1-1`, :math:`k_2 = 0,1,\ldots,n_2-1` and :math:`k_3 = 0,1,\ldots,n_3-1`. (Note the scale factor of :math:`\frac{1}{\sqrt{n_1n_2n_3}}` in this definition.)
The transformed values :math:`\hat{z}_{{k_1k_2k_3}}` are complex.
Because of conjugate symmetry (i.e., :math:`\hat{z}_{{k_1k_2k_3}}` is the complex conjugate of :math:`\hat{z}_{{\left(n_1-k_1\right)\left(n_2-k_2\right)\left(n_3-k_3\right)}}`), only slightly more than half of the Fourier coefficients need to be stored in the output.
A call of ``fft_real_3d`` followed by a call of :meth:`fft_hermitian_3d` will restore the original data.
This function calls :meth:`fft_realherm_1d_multi_col` and :meth:`fft_complex_1d_multi_row` to perform multiple one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm in Brigham (1974) and Temperton (1983).
.. _c06py-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_hermitian_3d(n1, n2, n3, y):
r"""
``fft_hermitian_3d`` computes the three-dimensional inverse discrete Fourier transform of a trivariate Hermitian sequence of complex data values.
.. _c06pz-py2-py-doc:
For full information please refer to the NAG Library document for c06pz
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06pzf.html
.. _c06pz-py2-py-parameters:
**Parameters**
**n1** : int
:math:`n_1`, the first dimension of the transform.
**n2** : int
:math:`n_2`, the second dimension of the transform.
**n3** : int
:math:`n_3`, the third dimension of the transform.
**y** : complex, array-like, shape :math:`\left(\left(\mathrm{n1}/2+1\right)\times \mathrm{n2}\times \mathrm{n3}\right)`
The Hermitian sequence of complex input dataset :math:`z`, where :math:`z_{{\textit{j}_1\textit{j}_2\textit{j}_3}}` is stored in :math:`\mathrm{y}[ \textit{j}_3 \times \left(n_1/2+1\right) n_2 + \textit{j}_2 \times \left(n_1/2+1\right) + \textit{j}_1 +1 -1]`, for :math:`\textit{j}_3 = 0,1,\ldots,n_3-1`, for :math:`\textit{j}_2 = 0,1,\ldots,n_2-1`, for :math:`\textit{j}_1 = 0,1,\ldots,n_1/2`.
**Returns**
**x** : float, ndarray, shape :math:`\left(\mathrm{n1}\times \mathrm{n2}\times \mathrm{n3}\right)`
The real output dataset :math:`\hat{x}`, where :math:`\hat{x}_{{\textit{k}_1\textit{k}_2\textit{k}_3}}` is stored in :math:`\mathrm{x}[ \textit{k}_3 \times n_1 n_2 + \textit{k}_2 \times n_1 + \textit{k}_1 +1 -1]`, for :math:`\textit{k}_3 = 0,1,\ldots,n_3-1`, for :math:`\textit{k}_2 = 0,1,\ldots,n_2-1`, for :math:`\textit{k}_1 = 0,1,\ldots,n_1-1`.
.. _c06pz-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{n1} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n1} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n2} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n2} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{n3} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n3} \geq 1`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06pz-py2-py-notes:
**Notes**
``fft_hermitian_3d`` computes the three-dimensional inverse discrete Fourier transform of a trivariate Hermitian sequence of complex data values :math:`z_{{j_1j_2j_3}}`, for :math:`\textit{j}_3 = 0,1,\ldots,n_3-1`, for :math:`\textit{j}_2 = 0,1,\ldots,n_2-1`, for :math:`\textit{j}_1 = 0,1,\ldots,n_1-1`.
The discrete Fourier transform is here defined by
.. math::
\hat{x}_{{k_1k_2k_3}} = \frac{1}{{\sqrt{n_1n_2n_3}}}\sum_{0}^{n_1-1}{\sum_{0}^{n_2-1}{\sum_{0}^{n_3-1}{z_{{j_1j_2j_3}}}}}\times \mathrm{exp}\left(2πi\left(\frac{{j_1k_1}}{{n_1}}+\frac{{j_2k_2}}{{n_2}}+\frac{{j_3k_3}}{{n_3}}\right)\right)\text{,}
where :math:`k_1 = 0,1,\ldots,n_1-1`, :math:`k_2 = 0,1,\ldots,n_2-1` and :math:`k_3 = 0,1,\ldots,n_3-1`. (Note the scale factor of :math:`\frac{1}{\sqrt{n_1n_2n_3}}` in this definition.)
Because the input data satisfies conjugate symmetry (i.e., :math:`z_{{j_1j_2j_3}}` is the complex conjugate of :math:`z_{{\left(n_1-j_1\right)\left(n_2-j_2\right)\left(n_3-j_3\right)}}`), the transformed values :math:`\hat{x}_{{k_1k_2k_3}}` are real.
A call of :meth:`fft_real_3d` followed by a call of ``fft_hermitian_3d`` will restore the original data.
This function calls :meth:`fft_realherm_1d_multi_col` and :meth:`fft_complex_1d_multi_row` to perform multiple one-dimensional discrete Fourier transforms by the fast Fourier transform (FFT) algorithm in Brigham (1974) and Temperton (1983).
.. _c06pz-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_real_sine_simple(m, n, x):
r"""
``fft_real_sine_simple`` computes the discrete Fourier sine transforms of :math:`m` sequences of real data values.
.. _c06ra-py2-py-doc:
For full information please refer to the NAG Library document for c06ra
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06raf.html
.. _c06ra-py2-py-parameters:
**Parameters**
**m** : int
:math:`m`, the number of sequences to be transformed.
**n** : int
One more than the number of real values in each sequence, i.e., the number of values in each sequence is :math:`n-1`.
**x** : float, array-like, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+2\right)\right)`
The data must be stored in :math:`\mathrm{x}` as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {1:\mathrm{n}+2}\right)`; each of the :math:`m` sequences is stored in a **row** of the array. In other words, if the :math:`n-1` data values of the :math:`\textit{p}`\ th sequence to be transformed are denoted by :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,{n-1}`, the first :math:`m\left(n-1\right)` elements of the array :math:`\mathrm{x}` must contain the values
.. math::
x_1^1,x_1^2,\ldots,x_1^m,x_2^1,x_2^2,\ldots,x_2^m,\ldots,x_{{n-1}}^1,x_{{n-1}}^2,\ldots,x_{{n-1}}^m\text{.}
The :math:`n`\ th to :math:`\left(n+2\right)`\ th elements of each row :math:`x_n^{\textit{p}},\ldots,x_{{n+2}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, are required as workspace. These :math:`3m` elements may contain arbitrary values as they are set to zero by the function.
**Returns**
**x** : float, ndarray, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+2\right)\right)`
The :math:`m` Fourier sine transforms stored as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {1:\mathrm{n}+2}\right)`. Each of the :math:`m` transforms is stored in a **row** of the array, overwriting the corresponding original sequence. If the :math:`\left(n-1\right)` components of the :math:`p`\ th Fourier sine transform are denoted by :math:`\hat{x}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 1,2,\ldots,{n-1}`, the :math:`m\left(n+2\right)` elements of the array :math:`\mathrm{x}` contain the values
.. math::
\hat{x}_1^1,\hat{x}_1^2,\ldots,\hat{x}_1^m,\hat{x}_2^1,\hat{x}_2^2,\ldots,\hat{x}_2^m,\ldots,\hat{x}_{{n-1}}^1,\hat{x}_{{n-1}}^2,\ldots,\hat{x}_{{n-1}}^m,0,0,\ldots,0\text{ }\left(3m\text{ times}\right)\text{.}
.. _c06ra-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06ra-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given :math:`m` sequences of :math:`{n-1}` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,{n-1}`, ``fft_real_sine_simple`` simultaneously calculates the Fourier sine transforms of all the sequences defined by
.. math::
\hat{x}_k^p = \sqrt{\frac{2}{n}}\sum_{{j = 1}}^{{n-1}}x_j^p\times \sin\left({jk}\frac{\pi }{n}\right)\text{, }\quad k = 1,2,\ldots,{n-1}\text{ and }p = 1,2,\ldots,m\text{.}
(Note the scale factor :math:`\sqrt{\frac{2}{n}}` in this definition.)
Since the Fourier sine transform defined above is its own inverse, two consecutive calls of this function will restore the original data.
The transform calculated by this function can be used to solve Poisson's equation when the solution is specified at both left and right boundaries (see Swarztrauber (1977)).
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, described in Temperton (1983), together with pre - and post-processing stages described in Swarztrauber (1982).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06ra-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Swarztrauber, P N, 1977, `The methods of cyclic reduction, Fourier analysis and the FACR algorithm for the discrete solution of Poisson's equation on a rectangle`, SIAM Rev. (19(3)), 490--501
Swarztrauber, P N, 1982, `Vectorizing the FFT's`, Parallel Computation, (ed G Rodrique), 51--83, Academic Press
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_real_cosine_simple(m, n, x):
r"""
``fft_real_cosine_simple`` computes the discrete Fourier cosine transforms of :math:`m` sequences of real data values.
.. _c06rb-py2-py-doc:
For full information please refer to the NAG Library document for c06rb
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06rbf.html
.. _c06rb-py2-py-parameters:
**Parameters**
**m** : int
:math:`m`, the number of sequences to be transformed.
**n** : int
One less than the number of real values in each sequence, i.e., the number of values in each sequence is :math:`n+1`.
**x** : float, array-like, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+3\right)\right)`
The data must be stored in :math:`\mathrm{x}` as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {0:\mathrm{n}+2}\right)`; each of the :math:`m` sequences is stored in a **row** of the array. In other words, if the :math:`\left(n+1\right)` data values of the :math:`\textit{p}`\ th sequence to be transformed are denoted by :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,n`, the first :math:`m\left(n+1\right)` elements of the array :math:`\mathrm{x}` must contain the values
.. math::
x_0^1,x_0^2,\ldots,x_0^m,x_1^1,x_1^2,\ldots,x_1^m,\ldots,x_n^1,x_n^2,\ldots,x_n^m\text{.}
The :math:`\left(n+2\right)`\ th and :math:`\left(n+3\right)`\ th elements of each row :math:`x_{{n+1}}^{\textit{p}},x_{{n+2}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, are required as workspace. These :math:`2m` elements may contain arbitrary values as they are set to zero by the function.
**Returns**
**x** : float, ndarray, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+3\right)\right)`
The :math:`m` Fourier cosine transforms stored as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {0:\mathrm{n}+2}\right)`. Each of the :math:`m` transforms is stored in a **row** of the array, overwriting the corresponding original data. If the :math:`\left(n+1\right)` components of the :math:`\textit{p}`\ th Fourier cosine transform are denoted by :math:`\hat{x}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 0,1,\ldots,n`, the :math:`m\left(n+3\right)` elements of the array :math:`\mathrm{x}` contain the values
.. math::
\hat{x}_0^1,\hat{x}_0^2,\ldots,\hat{x}_0^m,\hat{x}_1^1,\hat{x}_1^2,\ldots,\hat{x}_1^m,\ldots,\hat{x}_n^1,\hat{x}_n^2,\ldots,\hat{x}_n^m,0,0,\ldots,0\text{ }\left(2m\text{ times}\right)\text{.}
.. _c06rb-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06rb-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given :math:`m` sequences of :math:`n+1` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,n`, ``fft_real_cosine_simple`` simultaneously calculates the Fourier cosine transforms of all the sequences defined by
.. math::
\hat{x}_k^p = \sqrt{\frac{2}{n}}\left(\frac{1}{2}x_0^p+\sum_{{j = 1}}^{{n-1}}x_j^p\times \cos\left({jk}\frac{\pi }{n}\right)+\frac{1}{2}\left(-1\right)^kx_n^p\right)\text{, }\quad k = 0,1,\ldots,n\text{ and }p = 1,2,\ldots,m\text{.}
(Note the scale factor :math:`\sqrt{\frac{2}{n}}` in this definition.)
Since the Fourier cosine transform is its own inverse, two consecutive calls of ``fft_real_cosine_simple`` will restore the original data.
The transform calculated by this function can be used to solve Poisson's equation when the derivative of the solution is specified at both left and right boundaries (see Swarztrauber (1977)).
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, described in Temperton (1983), together with pre - and post-processing stages described in Swarztrauber (1982).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06rb-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Swarztrauber, P N, 1977, `The methods of cyclic reduction, Fourier analysis and the FACR algorithm for the discrete solution of Poisson's equation on a rectangle`, SIAM Rev. (19(3)), 490--501
Swarztrauber, P N, 1982, `Vectorizing the FFT's`, Parallel Computation, (ed G Rodrique), 51--83, Academic Press
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_real_qtrsine_simple(direct, m, n, x):
r"""
``fft_real_qtrsine_simple`` computes the discrete quarter-wave Fourier sine transforms of :math:`m` sequences of real data values.
.. _c06rc-py2-py-doc:
For full information please refer to the NAG Library document for c06rc
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06rcf.html
.. _c06rc-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06rc-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**m** : int
:math:`m`, the number of sequences to be transformed.
**n** : int
:math:`n`, the number of real values in each sequence.
**x** : float, array-like, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+2\right)\right)`
The data must be stored in :math:`\mathrm{x}` as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {1:\mathrm{n}+2}\right)`; each of the :math:`m` sequences is stored in a **row** of the array. In other words, if the data values of the :math:`\textit{p}`\ th sequence to be transformed are denoted by :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,n`, the first :math:`mn` elements of the array :math:`\mathrm{x}` must contain the values
.. math::
x_1^1,x_1^2,\ldots,x_1^m,x_2^1,x_2^2,\ldots,x_2^m,\ldots,x_n^1,x_n^2,\ldots,x_n^m\text{.}
The :math:`\left(n+1\right)`\ th and :math:`\left(n+2\right)`\ th elements of each row :math:`x_{{n+1}}^{\textit{p}},x_{{n+2}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, are required as workspace. These :math:`2m` elements may contain arbitrary values as they are set to zero by the function.
**Returns**
**x** : float, ndarray, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+2\right)\right)`
The :math:`m` quarter-wave sine transforms stored as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {1:\mathrm{n}+2}\right)`. Each of the :math:`m` transforms is stored in a **row** of the array, overwriting the corresponding original sequence. If the :math:`n` components of the :math:`\textit{p}`\ th quarter-wave sine transform are denoted by :math:`\hat{x}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 1,2,\ldots,n`, the :math:`m\left(n+2\right)` elements of the array :math:`\mathrm{x}` contain the values
.. math::
\hat{x}_1^1,\hat{x}_1^2,\ldots,\hat{x}_1^m,\hat{x}_2^1,\hat{x}_2^2,\ldots,\hat{x}_2^m,\ldots,\hat{x}_n^1,\hat{x}_n^2,\ldots,\hat{x}_n^m,0,0,\ldots,0\text{ }\left(2m\text{ times}\right)\text{.}
.. _c06rc-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06rc-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given :math:`m` sequences of :math:`n` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,n`, ``fft_real_qtrsine_simple`` simultaneously calculates the quarter-wave Fourier sine transforms of all the sequences defined by
.. math::
\hat{x}_k^p = \frac{1}{\sqrt{n}}\left(\sum_{{j = 1}}^{{n-1}}x_j^p\times \sin\left(j\left(2k-1\right)\frac{\pi }{{2n}}\right)+\frac{1}{2}\left(-1\right)^{{k-1}}x_n^p\right)\text{, if }\mathrm{direct} = \texttt{'F'}\text{,}
or its inverse
.. math::
x_k^p = \frac{2}{\sqrt{n}}\sum_{{j = 1}}^n\hat{x}_j^p\times \sin\left(\left(2j-1\right)k\frac{\pi }{{2n}}\right)\text{, if }\mathrm{direct} = \texttt{'B'}\text{,}
where :math:`k = 1,2,\ldots,n` and :math:`p = 1,2,\ldots,m`.
(Note the scale factor :math:`\frac{1}{\sqrt{n}}` in this definition.)
A call of ``fft_real_qtrsine_simple`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
The transform calculated by this function can be used to solve Poisson's equation when the solution is specified at the left boundary, and the derivative of the solution is specified at the right boundary (see Swarztrauber (1977)).
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, described in Temperton (1983), together with pre - and post-processing stages described in Swarztrauber (1982).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06rc-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Swarztrauber, P N, 1977, `The methods of cyclic reduction, Fourier analysis and the FACR algorithm for the discrete solution of Poisson's equation on a rectangle`, SIAM Rev. (19(3)), 490--501
Swarztrauber, P N, 1982, `Vectorizing the FFT's`, Parallel Computation, (ed G Rodrique), 51--83, Academic Press
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_real_qtrcosine_simple(direct, m, n, x):
r"""
``fft_real_qtrcosine_simple`` computes the discrete quarter-wave Fourier cosine transforms of :math:`m` sequences of real data values.
.. _c06rd-py2-py-doc:
For full information please refer to the NAG Library document for c06rd
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06rdf.html
.. _c06rd-py2-py-parameters:
**Parameters**
**direct** : str, length 1
If the forward transform as defined in :ref:`Notes <c06rd-py2-py-notes>` is to be computed, :math:`\mathrm{direct}` must be set equal to 'F'.
If the backward transform is to be computed, :math:`\mathrm{direct}` must be set equal to 'B'.
**m** : int
:math:`m`, the number of sequences to be transformed.
**n** : int
:math:`n`, the number of real values in each sequence.
**x** : float, array-like, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+2\right)\right)`
The data must be stored in :math:`\mathrm{x}` as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {0:\mathrm{n}+1}\right)`; each of the :math:`m` sequences is stored in a **row** of the array. In other words, if the data values of the :math:`\textit{p}`\ th sequence to be transformed are denoted by :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, the first :math:`mn` elements of the array :math:`\mathrm{x}` must contain the values
.. math::
x_0^1,x_0^2,\ldots,x_0^m,x_1^1,x_1^2,\ldots,x_1^m,\ldots,x_{{n-1}}^1,x_{{n-1}}^2,\ldots,x_{{n-1}}^m\text{.}
The :math:`\left(n+1\right)`\ th and :math:`\left(n+2\right)`\ th elements of each row :math:`x_n^{\textit{p}},x_{{n+1}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, are required as workspace. These :math:`2m` elements may contain arbitrary values as they are set to zero by the function.
**Returns**
**x** : float, ndarray, shape :math:`\left(\mathrm{m}\times \left(\mathrm{n}+2\right)\right)`
The :math:`m` quarter-wave cosine transforms stored as if in a two-dimensional array of dimension :math:`\left({1:\mathrm{m}}, {0:\mathrm{n}+1}\right)`. Each of the :math:`m` transforms is stored in a **row** of the array, overwriting the corresponding original sequence. If the :math:`n` components of the :math:`\textit{p}`\ th quarter-wave cosine transform are denoted by :math:`\hat{x}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 0,1,\ldots,{n-1}`, the :math:`m\left(n+2\right)` elements of the array :math:`\mathrm{x}` contain the values
.. math::
\hat{x}_0^1,\hat{x}_0^2,\ldots,\hat{x}_0^m,\hat{x}_1^1,\hat{x}_1^2,\ldots,\hat{x}_1^m,\ldots,\hat{x}_{{n-1}}^1,\hat{x}_{{n-1}}^2,\ldots,\hat{x}_{{n-1}}^m,0,0,\ldots,0\text{ }\left(2m\text{ times}\right)\text{.}
.. _c06rd-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{m} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{m} \geq 1`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{n} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{n} \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{direct} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{direct} = \texttt{'F'}` or :math:`\texttt{'B'}`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06rd-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
Given :math:`m` sequences of :math:`n` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_real_qtrcosine_simple`` simultaneously calculates the quarter-wave Fourier cosine transforms of all the sequences defined by
.. math::
\hat{x}_k^p = \frac{1}{\sqrt{n}}\left(\frac{1}{2}x_0^p+\sum_{{j = 1}}^{{n-1}}x_j^p\times \cos\left(j\left(2k-1\right)\frac{\pi }{{2n}}\right)\right)\text{, if }\mathrm{direct} = \texttt{'F'}\text{,}
or its inverse
.. math::
x_k^p = \frac{2}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}\hat{x}_j^p\times \cos\left(\left(2j-1\right)k\frac{\pi }{{2n}}\right)\text{, if }\mathrm{direct} = \texttt{'B'}\text{,}
where :math:`k = 0,1,\ldots,{n-1}` and :math:`p = 1,2,\ldots,m`.
(Note the scale factor :math:`\frac{1}{\sqrt{n}}` in this definition.)
A call of ``fft_real_qtrcosine_simple`` with :math:`\mathrm{direct} = \texttt{'F'}` followed by a call with :math:`\mathrm{direct} = \texttt{'B'}` will restore the original data.
The transform calculated by this function can be used to solve Poisson's equation when the derivative of the solution is specified at the left boundary, and the solution is specified at the right boundary (see Swarztrauber (1977)).
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, described in Temperton (1983), together with pre - and post-processing stages described in Swarztrauber (1982).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06rd-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Swarztrauber, P N, 1977, `The methods of cyclic reduction, Fourier analysis and the FACR algorithm for the discrete solution of Poisson's equation on a rectangle`, SIAM Rev. (19(3)), 490--501
Swarztrauber, P N, 1982, `Vectorizing the FFT's`, Parallel Computation, (ed G Rodrique), 51--83, Academic Press
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_sine(x):
r"""
``fft_sine`` computes the discrete Fourier sine transforms of :math:`m` sequences of real data values.
The elements of each sequence and its transform are stored contiguously.
.. _c06re-py2-py-doc:
For full information please refer to the NAG Library document for c06re
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06ref.html
.. _c06re-py2-py-parameters:
**Parameters**
**x** : float, array-like, shape :math:`\left(n-1, m\right)`
The data values of the :math:`\textit{p}`\ th sequence to be transformed, denoted by :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,{n-1}`, must be stored in :math:`\mathrm{x}[j-1,p-1]`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n-1, m\right)`
The :math:`\left(n-1\right)` components of the :math:`\textit{p}`\ th Fourier sine transform, denoted by :math:`\hat{x}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 1,2,\ldots,{n-1}`, are stored in :math:`\mathrm{x}[k-1,p-1]`, overwriting the corresponding original values.
.. _c06re-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`m \geq 1`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`3`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06re-py2-py-notes:
**Notes**
Given :math:`m` sequences of :math:`n-1` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,{n-1}`, ``fft_sine`` simultaneously calculates the Fourier sine transforms of all the sequences defined by
.. math::
\hat{x}_k^p = \sqrt{\frac{2}{n}}\sum_{{j = 1}}^{{n-1}}x_j^p\times \sin\left({jk}\frac{\pi }{n}\right)\text{, }\quad k = 1,2,\ldots,{n-1}\text{ and }p = 1,2,\ldots,m\text{.}
(Note the scale factor :math:`\sqrt{\frac{2}{n}}` in this definition.)
This transform is also known as type-I DST.
Since the Fourier sine transform defined above is its own inverse, two consecutive calls of ``fft_sine`` will restore the original data.
The transform calculated by this function can be used to solve Poisson's equation when the solution is specified at both left and right boundaries (see Swarztrauber (1977)).
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, described in Temperton (1983), together with pre - and post-processing stages described in Swarztrauber (1982).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06re-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Swarztrauber, P N, 1977, `The methods of cyclic reduction, Fourier analysis and the FACR algorithm for the discrete solution of Poisson's equation on a rectangle`, SIAM Rev. (19(3)), 490--501
Swarztrauber, P N, 1982, `Vectorizing the FFT's`, Parallel Computation, (ed G Rodrique), 51--83, Academic Press
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_cosine(x):
r"""
``fft_cosine`` computes the discrete Fourier cosine transforms of :math:`m` sequences of real data values.
The elements of each sequence and its transform are stored contiguously.
.. _c06rf-py2-py-doc:
For full information please refer to the NAG Library document for c06rf
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06rff.html
.. _c06rf-py2-py-parameters:
**Parameters**
**x** : float, array-like, shape :math:`\left(n+1, m\right)`
The data values of the :math:`\textit{p}`\ th sequence to be transformed, denoted by :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,n`, must be stored in :math:`\mathrm{x}[j,p-1]`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n+1, m\right)`
The :math:`\left(n+1\right)` components of the :math:`\textit{p}`\ th Fourier cosine transform, denoted by :math:`\hat{x}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 0,1,\ldots,n`, are stored in :math:`\mathrm{x}[k,p-1]`, overwriting the corresponding original values.
.. _c06rf-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`m \geq 1`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`3`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06rf-py2-py-notes:
**Notes**
Given :math:`m` sequences of :math:`n+1` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,n`, ``fft_cosine`` simultaneously calculates the Fourier cosine transforms of all the sequences defined by
.. math::
\hat{x}_k^p = \sqrt{\frac{2}{n}}\left(\frac{1}{2}x_0^p+\sum_{{j = 1}}^{{n-1}}x_j^p\times \cos\left({jk}\frac{\pi }{n}\right)+\frac{1}{2}\left(-1\right)^kx_n^p\right)\text{, }\quad k = 0,1,\ldots,n\text{ and }p = 1,2,\ldots,m\text{.}
(Note the scale factor :math:`\sqrt{\frac{2}{n}}` in this definition.)
This transform is also known as type-I DCT.
Since the Fourier cosine transform defined above is its own inverse, two consecutive calls of ``fft_cosine`` will restore the original data.
The transform calculated by this function can be used to solve Poisson's equation when the derivative of the solution is specified at both left and right boundaries (see Swarztrauber (1977)).
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, described in Temperton (1983), together with pre - and post-processing stages described in Swarztrauber (1982).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06rf-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Swarztrauber, P N, 1977, `The methods of cyclic reduction, Fourier analysis and the FACR algorithm for the discrete solution of Poisson's equation on a rectangle`, SIAM Rev. (19(3)), 490--501
Swarztrauber, P N, 1982, `Vectorizing the FFT's`, Parallel Computation, (ed G Rodrique), 51--83, Academic Press
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_qtrsine(idir, x):
r"""
``fft_qtrsine`` computes the discrete quarter-wave Fourier sine transforms of :math:`m` sequences of real data values.
The elements of each sequence and its transform are stored contiguously.
.. _c06rg-py2-py-doc:
For full information please refer to the NAG Library document for c06rg
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06rgf.html
.. _c06rg-py2-py-parameters:
**Parameters**
**idir** : int
Indicates the transform, as defined in :ref:`Notes <c06rg-py2-py-notes>`, to be computed.
:math:`\mathrm{idir} = 1`
Forward transform.
:math:`\mathrm{idir} = -1`
Inverse transform.
**x** : float, array-like, shape :math:`\left(n, m\right)`
The data values of the :math:`\textit{p}`\ th sequence to be transformed, denoted by :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,n`, must be stored in :math:`\mathrm{x}[j-1,p-1]`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n, m\right)`
The :math:`n` components of the :math:`\textit{p}`\ th quarter-wave sine transform, denoted by :math:`\hat{x}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 1,2,\ldots,n`, are stored in :math:`\mathrm{x}[k-1,p-1]`, overwriting the corresponding original values.
.. _c06rg-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`m \geq 1`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{idir} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{idir} = -1` or :math:`1`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06rg-py2-py-notes:
**Notes**
Given :math:`m` sequences of :math:`n` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 1,2,\ldots,n`, ``fft_qtrsine`` simultaneously calculates the quarter-wave Fourier sine transforms of all the sequences defined by
.. math::
\hat{x}_k^p = \frac{1}{\sqrt{n}}\left(\sum_{{j = 1}}^{{n-1}}x_j^p\times \sin\left(j\left(2k-1\right)\frac{\pi }{{2n}}\right)+\frac{1}{2}\left(-1\right)^{{k-1}}x_n^p\right)\text{, if }\mathrm{idir} = 1\text{,}
or its inverse
.. math::
x_k^p = \frac{2}{\sqrt{n}}\sum_{{j = 1}}^n\hat{x}_j^p\times \sin\left(\left(2j-1\right)k\frac{\pi }{{2n}}\right)\text{, if }\mathrm{idir} = -1\text{,}
where :math:`k = 1,2,\ldots,n` and :math:`p = 1,2,\ldots,m`.
(Note the scale factor :math:`\frac{1}{\sqrt{n}}` in this definition.)
A call of ``fft_qtrsine`` with :math:`\mathrm{idir} = 1` followed by a call with :math:`\mathrm{idir} = -1` will restore the original data.
The two transforms are also known as type-III DST and type-II DST, respectively.
The transform calculated by this function can be used to solve Poisson's equation when the solution is specified at the left boundary, and the derivative of the solution is specified at the right boundary (see Swarztrauber (1977)).
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, described in Temperton (1983), together with pre - and post-processing stages described in Swarztrauber (1982).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06rg-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Swarztrauber, P N, 1977, `The methods of cyclic reduction, Fourier analysis and the FACR algorithm for the discrete solution of Poisson's equation on a rectangle`, SIAM Rev. (19(3)), 490--501
Swarztrauber, P N, 1982, `Vectorizing the FFT's`, Parallel Computation, (ed G Rodrique), 51--83, Academic Press
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fft_qtrcosine(idir, x):
r"""
``fft_qtrcosine`` computes the discrete quarter-wave Fourier cosine transforms of :math:`m` sequences of real data values.
The elements of each sequence and its transform are stored contiguously.
.. _c06rh-py2-py-doc:
For full information please refer to the NAG Library document for c06rh
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06rhf.html
.. _c06rh-py2-py-parameters:
**Parameters**
**idir** : int
Indicates the transform, as defined in :ref:`Notes <c06rh-py2-py-notes>`, to be computed.
:math:`\mathrm{idir} = 1`
Forward transform.
:math:`\mathrm{idir} = -1`
Inverse transform.
**x** : float, array-like, shape :math:`\left(n, m\right)`
The data values of the :math:`\textit{p}`\ th sequence to be transformed, denoted by :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,n-1`, must be stored in :math:`\mathrm{x}[j,p-1]`.
**Returns**
**x** : float, ndarray, shape :math:`\left(n, m\right)`
The :math:`n` components of the :math:`\textit{p}`\ th quarter-wave cosine transform, denoted by :math:`\hat{x}_{\textit{k}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{k} = 0,1,\ldots,n-1`, are stored in :math:`\mathrm{x}[k,p-1]`, overwriting the corresponding original values.
.. _c06rh-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`m \geq 1`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{idir} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{idir} = -1` or :math:`1`.
(`errno` :math:`4`)
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact `NAG <https://www.nag.com>`__ for assistance.
.. _c06rh-py2-py-notes:
**Notes**
Given :math:`m` sequences of :math:`n` real data values :math:`x_{\textit{j}}^{\textit{p}}`, for :math:`\textit{p} = 1,2,\ldots,m`, for :math:`\textit{j} = 0,1,\ldots,{n-1}`, ``fft_qtrcosine`` simultaneously calculates the quarter-wave Fourier cosine transforms of all the sequences defined by
.. math::
\hat{x}_k^p = \frac{1}{\sqrt{n}}\left(\frac{1}{2}x_0^p+\sum_{{j = 1}}^{{n-1}}x_j^p\times \cos\left(j\left(2k+1\right)\frac{\pi }{{2n}}\right)\right)\text{, if }\mathrm{idir} = 1\text{,}
or its inverse
.. math::
x_k^p = \frac{2}{\sqrt{n}}\sum_{{j = 0}}^{{n-1}}\hat{x}_j^p\times \cos\left(\left(2j+1\right)k\frac{\pi }{{2n}}\right)\text{, if }\mathrm{idir} = -1\text{,}
where :math:`k = 0,1,\ldots,{n-1}` and :math:`p = 1,2,\ldots,m`.
(Note the scale factor :math:`\frac{1}{\sqrt{n}}` in this definition.)
A call of ``fft_qtrcosine`` with :math:`\mathrm{idir} = 1` followed by a call with :math:`\mathrm{idir} = -1` will restore the original data.
The two transforms are also known as type-III DCT and type-II DCT, respectively.
The transform calculated by this function can be used to solve Poisson's equation when the derivative of the solution is specified at the left boundary, and the solution is specified at the right boundary (see Swarztrauber (1977)).
The function uses a variant of the fast Fourier transform (FFT) algorithm (see Brigham (1974)) known as the Stockham self-sorting algorithm, described in Temperton (1983), together with pre - and post-processing stages described in Swarztrauber (1982).
Special coding is provided for the factors :math:`2`, :math:`3`, :math:`4` and :math:`5`.
.. _c06rh-py2-py-references:
**References**
Brigham, E O, 1974, `The Fast Fourier Transform`, Prentice--Hall
Swarztrauber, P N, 1977, `The methods of cyclic reduction, Fourier analysis and the FACR algorithm for the discrete solution of Poisson's equation on a rectangle`, SIAM Rev. (19(3)), 490--501
Swarztrauber, P N, 1982, `Vectorizing the FFT's`, Parallel Computation, (ed G Rodrique), 51--83, Academic Press
Temperton, C, 1983, `Fast mixed-radix real Fourier transforms`, J. Comput. Phys. (52), 340--350
"""
raise NotImplementedError
[docs]def fast_gauss(srcs, trgs, q, p1, p2, k, hin, tol):
r"""
``fast_gauss`` calculates the multidimensional fast Gauss transform.
.. _c06sa-py2-py-doc:
For full information please refer to the NAG Library document for c06sa
https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06saf.html
.. _c06sa-py2-py-parameters:
**Parameters**
**srcs** : float, array-like, shape :math:`\left(d, n\right)`
:math:`x`, the locations of the Gaussian sources.
**trgs** : float, array-like, shape :math:`\left(d, m\right)`
:math:`y`, the locations of the target points at which the FGT will be evaluated.
**q** : float, array-like, shape :math:`\left(n\right)`
:math:`q`, the weights of the Gaussian sources.
**p1** : int
:math:`p_1`, the number of terms of the first Taylor series to be evaluated.
**p2** : int
:math:`p_2`, the number of terms of the second Taylor series to be evaluated.
**k** : int
:math:`k`, the number of clusters into which the source points will be aggregated.
**hin** : float, array-like, shape :math:`\left(\textit{lhin}\right)`
:math:`h`, the standard deviations of the Gaussian sources. If :math:`\textit{lhin} < n`, the array will be expanded automatically by repeating :math:`\mathrm{hin}` until it is of length :math:`\textit{n}`. See `the G01 Introduction <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/g01/g01intro.html#vectorizeddoc>`__ for further information.
**tol** : float
:math:`\epsilon`, the desired accuracy of the FGT approximation of the DGT. Determines the radius of the source clusters: the contribution of a source point to the FGT approximation at a target point is disregarded if the source is outside the corresponding cluster radius.
**Returns**
**p1** : int
:math:`\mathrm{p1}` is unchanged.
**p2** : int
:math:`\mathrm{p2}` is unchanged.
**k** : int
:math:`\mathrm{k}` is unchanged.
**v** : float, ndarray, shape :math:`\left(m\right)`
:math:`\hat{G}\left(y\right)`, the value of the FGT evaluated at :math:`y`.
**term** : float, ndarray, shape :math:`\left(m\right)`
:math:`\mathrm{term}[j-1]` contains the absolute value of the final Taylor series term that is largest, relative to the size of the sum of the corresponding series, across all clusters that contribute to the FGT at target point :math:`\mathrm{v}[j-1]`.
.. _c06sa-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`d = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`d > 0`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n > 0`.
(`errno` :math:`3`)
On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`m > 0`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{p1} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{p1} > 0`.
(`errno` :math:`5`)
On entry, :math:`\mathrm{p2} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{p2} > 0`.
(`errno` :math:`6`)
On entry, :math:`\mathrm{k} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`1\leq \mathrm{k}\leq n`.
(`errno` :math:`7`)
On entry, :math:`\mathrm{hin}[\langle\mathit{\boldsymbol{value}}\rangle] = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{hin}[\textit{i}] > 0.0`, for :math:`\textit{i} = 0,\ldots,\textit{lhin}-1`.
(`errno` :math:`8`)
On entry, :math:`\textit{lhin} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`1\leq \textit{lhin}\leq n`.
(`errno` :math:`9`)
On entry, :math:`\mathrm{tol} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{tol} > 0.0`.
**Warns**
**NagAlgorithmicWarning**
(`errno` :math:`10`)
On exit, :math:`\mathrm{p1} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{p2} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{k} = \langle\mathit{\boldsymbol{value}}\rangle`.
:math:`\mathrm{p1}`, :math:`\mathrm{p2}` or :math:`\mathrm{k}` may have been too small to calculate :math:`\mathrm{v}[m-1]` to the required accuracy :math:`\mathrm{tol}`.
.. _c06sa-py2-py-notes:
**Notes**
``fast_gauss`` calculates the :math:`d`-dimensional fast Gauss transform (FGT), :math:`\hat{G}\left(y\right)`, that approximates the discrete Gauss transform (DGT), :math:`G\left(y\right)`, evaluated at a set of target points :math:`y_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,m \in \mathbb{R}^d`.
The DGT is defined as:
.. math::
G\left(y_j\right) = \sum_{{i = 1}}^nq_ie^{{-\left\lVert y_j-x_i\right\rVert_2^2/h_i^2}}\text{,}\quad j = 1,\ldots,m
where :math:`x_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n \in \mathbb{R}^d`, are the Gaussian source points, :math:`q_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n \in \mathbb{R}^+`, are the source weights and :math:`h_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n \in \mathbb{R}^+`, are the source standard deviations (alternatively source scales or source bandwidths).
This function implements the improved FGT algorithm presented in Raykar and Duraiswami (2005).
The algorithm clusters the sources into :math:`k` distinct clusters and then computes two Taylor series approximations per cluster with :math:`p_1` and :math:`p_2` terms respectively.
You must provide :math:`p_1`, :math:`p_2` and :math:`k` when calling the function.
See `Accuracy <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/c06/c06saf.html#accuracy>`__ below for a further discussion on accuracy when choosing their values.
The input array :math:`\mathrm{hin}` of this function is designed to allow maximum flexibility in the supply of the standard deviation arguments by reusing, in a cyclic manner, elements of the array when it is less than :math:`n` elements long.
For example, if all Gaussian sources have the same standard deviation then it is only necessary to set :math:`\textit{lhin}` to :math:`1` and to provide the value of the standard deviation in :math:`\mathrm{hin}[0]`; the function will then automatically expand :math:`\mathrm{hin}` to be of length :math:`n`.
For further details please see `the G01 Introduction <https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/g01/g01intro.html#vectorizeddoc>`__.
.. _c06sa-py2-py-references:
**References**
Greengard, L and Strain, J, 1991, `The Fast Gauss Transform`, SIAM J. Sci. Statist. Comput. (12(1)), 79--94
Raykar, V C and Duraiswami, R, 2005, `Improved Fast Gauss Transform With Variable Source Scales`, University of Maryland Technical Report CS-TR-4727/UMIACS-TR-2005-34
"""
raise NotImplementedError
