# -*- coding: utf-8 -*-
r"""
Module Summary
--------------
Interfaces for the NAG Mark 30.0 `zeros` Chapter.
``zeros`` - Zeros of Polynomials
This module is concerned with computing the zeros of a polynomial with real or complex coefficients.
See Also
--------
``naginterfaces.library.examples.zeros`` :
This subpackage contains examples for the ``zeros`` module.
See also the :ref:`library_zeros_ex` subsection.
Functionality Index
-------------------
**All zeros of cubic**
complex coefficients: :meth:`cubic_complex`
real coefficients: :meth:`cubic_real`
**All zeros of polynomial**
complex coefficients
fast modified Laguerre's method: :meth:`poly_complex_fpml`
modified Laguerre's method: :meth:`poly_complex`
real coefficients
fast modified Laguerre's method: :meth:`poly_real_fpml`
modified Laguerre's method: :meth:`poly_real`
**All zeros of quadratic**
complex coefficients: :meth:`quadratic_complex`
real coefficients: :meth:`quadratic_real`
**All zeros of quartic**
complex coefficients: :meth:`quartic_complex`
real coefficients: :meth:`quartic_real`
For full information please refer to the NAG Library document
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02intro.html
"""
# NAG Copyright 2017-2024.
[docs]def poly_complex_fpml(a, itmax=30, polish=1):
r"""
``poly_complex_fpml`` finds all the roots of a complex polynomial equation, using a fourth-order convergent modification of Laguerre's method.
.. _c02aa-py2-py-doc:
For full information please refer to the NAG Library document for c02aa
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02aaf.html
.. _c02aa-py2-py-parameters:
**Parameters**
**a** : complex, array-like, shape :math:`\left(n+1\right)`
:math:`\mathrm{a}[\textit{i}]` must contain the coefficient of :math:`z^{{n-\textit{i}}}`, for :math:`\textit{i} = 0,1,\ldots,n`.
**itmax** : int, optional
The maximum number of iterations to be performed.
**polish** : int, optional
Specifies the polishing technique used to refine root approximations.
:math:`\mathrm{polish} = 0`
No polishing.
:math:`\mathrm{polish} = 1`
Single iteration of Newton's method.
:math:`\mathrm{polish} = 2`
Iterative refinement using the compensated Horner's method.
See also `Selecting a Polishing Process <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02aaf.html#fcomments_polish>`__.
**Returns**
**z** : complex, ndarray, shape :math:`\left(n\right)`
:math:`\mathrm{z}[\textit{j}-1]` holds the approximation of the root :math:`z_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,n`.
**berr** : float, ndarray, shape :math:`\left(n\right)`
:math:`\mathrm{berr}[\textit{j}-1]` holds the relative backward error, :math:`\eta \left(z_{\textit{j}}\right)`, for :math:`\textit{j} = 1,2,\ldots,n`.
**cond** : float, ndarray, shape :math:`\left(n\right)`
:math:`\mathrm{cond}[\textit{j}-1]` holds the condition number, :math:`\kappa \left(z_{\textit{j}}\right)`, for :math:`\textit{j} = 1,2,\ldots,n`.
**conv** : int, ndarray, shape :math:`\left(n\right)`
:math:`\mathrm{conv}[\textit{j}-1]` indicates the convergence status of the root approximation, :math:`z_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,n`.
:math:`\mathrm{conv}[j-1]\geq 0`
Successfully converged after :math:`\mathrm{conv}[j-1]` iterations.
:math:`\mathrm{conv}[j-1] = -1`
Failed to converge after :math:`\mathrm{itmax}` iterations.
:math:`\mathrm{conv}[j-1] = -2`
Overflow was encountered.
Note: if :math:`\mathrm{polish} = 2`, :math:`\mathrm{conv}` refers to convergence in the compensated polishing process.
.. _c02aa-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, the complex variable :math:`\mathrm{a}[0] = \left(0.0, 0.0\right)`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{itmax} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{itmax}\geq 1`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{polish} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{polish} = 0`, :math:`1` or :math:`2`.
(`errno` :math:`5`)
Convergence has failed for at least one root approximation.
(`errno` :math:`6`)
``poly_complex_fpml`` encountered overflow during at least one root approximation.
.. _c02aa-py2-py-notes:
**Notes**
``poly_complex_fpml`` attempts to find all the roots of the :math:`n`\ th degree complex polynomial equation
.. math::
P\left(z\right) = a_0z^n+a_1z^{{n-1}}+a_2z^{{n-2}} + \cdots +a_{{n-1}}z+a_n = 0\text{.}
The roots are located using a modified form of Laguerre's method, as implemented by Cameron (2018).
An implicit deflation strategy is employed, which allows for high accuracy even when solving high degree polynomials.
Linear (:math:`n = 1`) and quadratic (:math:`n = 2`) problems are obtained via the 'standard' closed formulae.
First, initial estimates of the roots are made using a method proposed by Bini (1996), which selects complex numbers along circles of suitable radii.
Updates to each root approximation :math:`z_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,n`, are then made using the iterative formula from Petković `et al.` (1997):
.. math::
\hat{z}_j = z_j-\frac{{n}}{{G_j\pm \sqrt{\left(n-1\right)\left(nH_j-G_j^2\right)}}}\text{,}
where
.. math::
G_j = \frac{{P^{\prime }\left(z_j\right)}}{{P\left(z_j\right)}}-\sum_{1\text{; }\textit{i}\neq j}^{n}{\frac{{1}}{{\left(z_j-z_{\textit{i}}\right)}}}\hspace{1em}\text{and}\hspace{1em}H_j = -\left(\frac{{P^{\prime }\left(z_j\right)}}{{P\left(z_j\right)}}\right)^{\prime }-\sum_{1\text{; }\textit{i}\neq j}^{n}{\frac{{1}}{{\left(z_j-z_{\textit{i}}\right)^2}}}\text{.}
The nearest :math:`\hat{z}_j` to the current approximation is used, by selecting the sign that maximizes the denominator of the correction term.
The subtraction of the sum terms when computing :math:`G_j` and :math:`H_j` determines the implicit deflation strategy of the modified Laguerre method.
The relative backward error of a root approximation :math:`z_j` is given by
.. math::
\eta \left(z_j\right) = \frac{{\left\lvert P\left(z_j\right)\right\rvert }}{{\alpha \left(z_j\right)}}\text{,}
where :math:`\alpha` is the perturbed polynomial
.. math::
\alpha \left(z\right) = \sum_{0}^{n}{\left(3.8\textit{k}+1\right)\left\lvert a_{{n-\textit{k}}}\right\rvert \left\lvert z\right\rvert^{\textit{k}}}\text{.}
A root approximation is deemed to have converged if :math:`\eta \left(z_j\right)\leq 2\text{ }\times` machine precision, at which point updates of that root cease.
If the stopping criterion holds, then the computed root is the exact root of a polynomial whose coefficients are no more perturbed than the floating-point computation of :math:`P\left(z_j\right)`.
The condition number of each root is also computed, as a measure of sensitivity to changes in the coefficients of the polynomial.
It is given by
.. math::
\kappa \left(z_j\right) = \frac{{\alpha \left(z_j\right)}}{{\left\lvert z_j\right\rvert \left\lvert P^{\prime }\left(z_j\right)\right\rvert }}\text{.}
Root approximations can be further refined with optional 'polishing' processes.
A simple polishing process is provided that carries out a single iteration of Newton's method, which proves quick yet often effective.
Alternatively, a compensated polishing process from Cameron and O'Neill (2019) can be applied.
This iterative method combines the implicit deflation of the modified Laguerre method, with the accuracy of evaluating polynomials and their derivatives using the compensated Horner's method from Graillat `et al.` (2005).
Compensated polishing yields approximations with a limiting accuracy as if computed in twice the working precision.
It is recommended that you read `Selecting a Polishing Process <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02aaf.html#fcomments_polish>`__ for advice on selecting an appropriate polishing process.
.. _c02aa-py2-py-references:
**References**
Bini, D A, 1996, `Numerical computation of polynomial zeros by means of Aberth's method`, Numerical Algorithms (13), 179--200, Springer US
Cameron, T R, 2018, `An effective implementation of a modified Laguerre method for the roots of a polynomial`, Numerical Algorithms, Springer US, https://doi.org/10.1007/s11075-018-0641-9
Cameron, T R and O'Neill, A, 2019, `On a compensated polishing technique for polynomial root solvers`, To be published,
Graillat, S, Louvet, N, and Langlois, P, 2005, `Compensated Horner scheme`, Technical Report, Université de Perpignan Via Domitia
Petković, M, Ilić, S, and Tričković, S, 1997, `A family of simultaneous zero-finding methods`, Computers & Mathematics with Applications (Volume 34) (10), 49--59, https://doi.org/10.1016/S0898-1221(97)00206-X
Wilkinson, J H, 1959, `The evaluation of the zeros of ill-conditioned polynomials. Part I`, Numerische Mathematik (Volume 1) (1), 150--166, Springer-Verlag
See Also
--------
:meth:`naginterfaces.library.examples.zeros.poly_complex_fpml_ex.main`
"""
raise NotImplementedError
[docs]def poly_real_fpml(a, itmax=30, polish=1):
r"""
``poly_real_fpml`` finds all the roots of a real polynomial equation, using a fourth-order convergent modification of Laguerre's method.
.. _c02ab-py2-py-doc:
For full information please refer to the NAG Library document for c02ab
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02abf.html
.. _c02ab-py2-py-parameters:
**Parameters**
**a** : float, array-like, shape :math:`\left(n+1\right)`
:math:`\mathrm{a}[\textit{i}]` must contain the coefficient of :math:`z^{{n-\textit{i}}}`, for :math:`\textit{i} = 0,1,\ldots,n`.
**itmax** : int, optional
The maximum number of iterations to be performed.
**polish** : int, optional
Specifies the polishing technique used to refine root approximations.
:math:`\mathrm{polish} = 0`
No polishing.
:math:`\mathrm{polish} = 1`
Single iteration of Newton's method.
:math:`\mathrm{polish} = 2`
Iterative refinement using the compensated Horner's method.
See also `Selecting a Polishing Process <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02abf.html#fcomments_polish>`__.
**Returns**
**z** : complex, ndarray, shape :math:`\left(n\right)`
:math:`\mathrm{z}[\textit{j}-1]` holds the approximation of the root :math:`z_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,n`.
**berr** : float, ndarray, shape :math:`\left(n\right)`
:math:`\mathrm{berr}[\textit{j}-1]` holds the relative backward error, :math:`\eta \left(z_{\textit{j}}\right)`, for :math:`\textit{j} = 1,2,\ldots,n`.
**cond** : float, ndarray, shape :math:`\left(n\right)`
:math:`\mathrm{cond}[\textit{j}-1]` holds the condition number, :math:`\kappa \left(z_{\textit{j}}\right)`, for :math:`\textit{j} = 1,2,\ldots,n`.
**conv** : int, ndarray, shape :math:`\left(n\right)`
:math:`\mathrm{conv}[\textit{j}-1]` indicates the convergence status of the root approximation, :math:`z_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,n`.
:math:`\mathrm{conv}[j-1]\geq 0`
Successfully converged after :math:`\mathrm{conv}[j-1]` iterations.
:math:`\mathrm{conv}[j-1] = -1`
Failed to converge after :math:`\mathrm{itmax}` iterations.
:math:`\mathrm{conv}[j-1] = -2`
Overflow was encountered.
Note: if :math:`\mathrm{polish} = 2`, :math:`\mathrm{conv}` refers to convergence in the compensated polishing process.
.. _c02ab-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, the real variable :math:`\mathrm{a}[0] = 0.0`.
(`errno` :math:`2`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{itmax} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{itmax}\geq 1`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{polish} = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`\mathrm{polish} = 0`, :math:`1` or :math:`2`.
(`errno` :math:`5`)
Convergence has failed for at least one root approximation.
(`errno` :math:`6`)
``poly_real_fpml`` encountered overflow during at least one root approximation.
.. _c02ab-py2-py-notes:
**Notes**
``poly_real_fpml`` attempts to find all the roots of the :math:`n`\ th degree real polynomial equation
.. math::
P\left(z\right) = a_0z^n+a_1z^{{n-1}}+a_2z^{{n-2}} + \cdots +a_{{n-1}}z+a_n = 0\text{.}
The roots are located using a modified form of Laguerre's method, as implemented by Cameron (2018).
``poly_real_fpml`` is a wrapper around the corresponding complex routine :meth:`poly_complex_fpml`.
The relative backward error of a root approximation :math:`z_j` is given by
.. math::
\eta \left(z_j\right) = \frac{{\left\lvert P\left(z_j\right)\right\rvert }}{{\alpha \left(z_j\right)}}\text{,}
where :math:`\alpha` is the perturbed polynomial
.. math::
\alpha \left(z\right) = \sum_{0}^{n}{\left(3.8\textit{k}+1\right)\left\lvert a_{{n-\textit{k}}}\right\rvert \left\lvert z\right\rvert^{\textit{k}}}\text{.}
A root approximation is deemed to have converged if :math:`\eta \left(z_j\right)\leq 2\text{ }\times` machine precision, at which point updates of that root cease.
If the stopping criterion holds, then the computed root is the exact root of a polynomial whose coefficients are no more perturbed than the floating-point computation of :math:`P\left(z_j\right)`.
The condition number of each root is also computed, as a measure of sensitivity to changes in the coefficients of the polynomial.
It is given by
.. math::
\kappa \left(z_j\right) = \frac{{\alpha \left(z_j\right)}}{{\left\lvert z_j\right\rvert \left\lvert P^{\prime }\left(z_j\right)\right\rvert }}\text{.}
Root approximations can be further refined with optional 'polishing' processes.
A simple polishing process is provided that carries out a single iteration of Newton's method, which proves quick yet often effective.
Alternatively, a compensated polishing process from Cameron and O'Neill (2019) can be applied.
This iterative method combines the implicit deflation of the modified Laguerre method, with the accuracy of evaluating polynomials and their derivatives using the compensated Horner's method from Graillat `et al.` (2005).
Compensated polishing yields approximations with a limiting accuracy as if computed in twice the working precision.
It is recommended that you read `Selecting a Polishing Process <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02abf.html#fcomments_polish>`__ for advice on selecting an appropriate polishing process.
.. _c02ab-py2-py-references:
**References**
Bini, D A, 1996, `Numerical computation of polynomial zeros by means of Aberth's method`, Numerical Algorithms (13), 179--200, Springer US
Cameron, T R, 2018, `An effective implementation of a modified Laguerre method for the roots of a polynomial`, Numerical Algorithms, Springer US, https://doi.org/10.1007/s11075-018-0641-9
Cameron, T R and O'Neill, A, 2019, `On a compensated polishing technique for polynomial root solvers`, To be published,
Graillat, S, Louvet, N, and Langlois, P, 2005, `Compensated Horner scheme`, Technical Report, Université de Perpignan Via Domitia
Petković, M, Ilić, S, and Tričković, S, 1997, `A family of simultaneous zero-finding methods`, Computers & Mathematics with Applications (Volume 34) (10), 49--59, https://doi.org/10.1016/S0898-1221(97)00206-X
Wilkinson, J H, 1959, `The evaluation of the zeros of ill-conditioned polynomials. Part I`, Numerische Mathematik (Volume 1) (1), 150--166, Springer-Verlag
"""
raise NotImplementedError
[docs]def poly_complex(a, scal=True):
r"""
``poly_complex`` finds all the roots of a complex polynomial equation, using a variant of Laguerre's method.
.. deprecated:: 27.1.0.0
``poly_complex`` is deprecated.
Please use :meth:`poly_complex_fpml` instead.
See also the :ref:`Replacement Calls <replace>` document.
.. _c02af-py2-py-doc:
For full information please refer to the NAG Library document for c02af
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02aff.html
.. _c02af-py2-py-parameters:
**Parameters**
**a** : complex, array-like, shape :math:`\left(n+1\right)`
:math:`\mathrm{Re}\left(\mathrm{a}[\textit{i}]\right)` and :math:`\mathrm{Im}\left(\mathrm{a}[\textit{i}]\right)` must contain the real and imaginary parts of :math:`a_{\textit{i}}` (i.e., the coefficient of :math:`z^{{n-\textit{i}}}`), for :math:`\textit{i} = 0,1,\ldots,n`.
**scal** : bool, optional
Indicates whether or not the polynomial is to be scaled. See `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02aff.html#fcomments>`__ for advice on when it may be preferable to set :math:`\mathrm{scal} = \mathbf{False}` and for a description of the scaling strategy.
**Returns**
**z** : complex, ndarray, shape :math:`\left(n\right)`
The roots.
.. _c02af-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{a}[0] = 0.0`.
Constraint: :math:`\mathrm{a}[0]\neq 0.0`.
(`errno` :math:`1`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`2`)
The iterative procedure has failed to converge. This error is very unlikely to occur. If it does, please contact `NAG <https://www.nag.com>`__ immediately, as some basic assumption for the arithmetic has been violated.
(`errno` :math:`3`)
``poly_complex`` cannot evaluate :math:`P\left(z\right)` near some of its zeros without underflow. If this message occurs please contact `NAG <https://www.nag.com>`__.
(`errno` :math:`3`)
``poly_complex`` cannot evaluate :math:`P\left(z\right)` near some of its zeros without overflow. If this message occurs please contact `NAG <https://www.nag.com>`__.
.. _c02af-py2-py-notes:
**Notes**
`In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.`
``poly_complex`` attempts to find all the roots of the :math:`n`\ th degree complex polynomial equation
.. math::
P\left(z\right) = a_0z^n+a_1z^{{n-1}}+a_2z^{{n-2}} + \cdots +a_{{n-1}}z+a_n = 0\text{.}
The roots are located using a modified form of Laguerre's method, originally proposed by Smith (1967).
The method of Laguerre (see Wilkinson (1965)) can be described by the iterative scheme
.. math::
L\left(z_k\right) = z_{{k+1}}-z_k = \frac{{{-n}P\left(z_k\right)}}{{P^{\prime }\left(z_k\right)\pm \sqrt{H\left(z_k\right)}}}\text{,}
where :math:`H\left(z_k\right) = \left(n-1\right)\left[\left(n-1\right){\left(P^{\prime }\left(z_k\right)\right)}^2-nP\left(z_k\right)P^{{\prime \prime }}\left(z_k\right)\right]` and :math:`z_0` is specified.
The sign in the denominator is chosen so that the modulus of the Laguerre step at :math:`z_k`, viz. :math:`\left\lvert L\left(z_k\right)\right\rvert`, is as small as possible.
The method can be shown to be cubically convergent for isolated roots (real or complex) and linearly convergent for multiple roots.
The function generates a sequence of iterates :math:`z_1,z_2,z_3,\ldots \text{}`, such that :math:`\left\lvert P\left(z_{{k+1}}\right)\right\rvert < \left\lvert P\left(z_k\right)\right\rvert` and ensures that :math:`z_{{k+1}}+L\left(z_{{k+1}}\right)` 'roughly' lies inside a circular region of radius :math:`\left\lvert F\right\rvert` about :math:`z_k` known to contain a zero of :math:`P\left(z\right)`; that is, :math:`\left\lvert L\left(z_{{k+1}}\right)\right\rvert \leq \left\lvert F\right\rvert`, where :math:`F` denotes the Fejér bound (see Marden (1966)) at the point :math:`z_k`.
Following Smith (1967), :math:`F` is taken to be :math:`\mathrm{min}\left(B, {1.445nR}\right)`, where :math:`B` is an upper bound for the magnitude of the smallest zero given by
.. math::
B = 1.0001\times \mathrm{min}\left({\sqrt{n}L\left(z_k\right)}, \left\lvert r_1\right\rvert, \left\lvert a_n/a_0\right\rvert^{{1/n}}\right)\text{,}
:math:`r_1` is the zero :math:`X` of smaller magnitude of the quadratic equation
.. math::
\frac{{P^{{\prime \prime }}\left(z_k\right)}}{{2n\left(n-1\right)}}X^2+\frac{{P^{\prime }\left(z_k\right)}}{n}X+\frac{1}{2}P\left(z_k\right) = 0
and the Cauchy lower bound :math:`R` for the smallest zero is computed (using Newton's Method) as the positive root of the polynomial equation
.. math::
\left\lvert a_0\right\rvert z^n+\left\lvert a_1\right\rvert z^{{n-1}}+\left\lvert a_2\right\rvert z^{{n-2}} + \cdots +\left\lvert a_{{n-1}}\right\rvert z-\left\lvert a_n\right\rvert = 0\text{.}
Starting from the origin, successive iterates are generated according to the rule :math:`z_{{k+1}} = z_k+L\left(z_k\right)`, for :math:`k = 1,2,3,\ldots \text{}`, and :math:`L\left(z_k\right)` is 'adjusted' so that :math:`\left\lvert P\left(z_{{k+1}}\right)\right\rvert < \left\lvert P\left(z_k\right)\right\rvert` and :math:`\left\lvert L\left(z_{{k+1}}\right)\right\rvert \leq \left\lvert F\right\rvert`.
The iterative procedure terminates if :math:`P\left(z_{{k+1}}\right)` is smaller in absolute value than the bound on the rounding error in :math:`P\left(z_{{k+1}}\right)` and the current iterate :math:`z_p = z_{{k+1}}` is taken to be a zero of :math:`P\left(z\right)`.
The deflated polynomial :math:`\tilde{P}\left(z\right) = P\left(z\right)/\left(z-z_p\right)` of degree :math:`n-1` is then formed, and the above procedure is repeated on the deflated polynomial until :math:`n < 3`, whereupon the remaining roots are obtained via the 'standard' closed formulae for a linear (:math:`n = 1`) or quadratic (:math:`n = 2`) equation.
Note that :meth:`quadratic_complex`, :meth:`cubic_complex` and :meth:`quartic_complex` can be used to obtain the roots of a quadratic, cubic (:math:`n = 3`) and quartic (:math:`n = 4`) polynomial, respectively.
.. _c02af-py2-py-references:
**References**
Marden, M, 1966, `Geometry of polynomials`, Mathematical Surveys (3), American Mathematical Society, Providence, RI
Smith, B T, 1967, `ZERPOL: a zero finding algorithm for polynomials using Laguerre's method`, Technical Report, Department of Computer Science, University of Toronto, Canada
Thompson, K W, 1991, `Error analysis for polynomial solvers`, Fortran Journal (Volume 3) (3), 10--13
Wilkinson, J H, 1965, `The Algebraic Eigenvalue Problem`, Oxford University Press, Oxford
"""
raise NotImplementedError
[docs]def poly_real(a, scal=True):
r"""
``poly_real`` finds all the roots of a real polynomial equation, using a variant of Laguerre's method.
.. deprecated:: 27.1.0.0
``poly_real`` is deprecated.
Please use :meth:`poly_real_fpml` instead.
See also the :ref:`Replacement Calls <replace>` document.
.. _c02ag-py2-py-doc:
For full information please refer to the NAG Library document for c02ag
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02agf.html
.. _c02ag-py2-py-parameters:
**Parameters**
**a** : float, array-like, shape :math:`\left(n+1\right)`
:math:`\mathrm{a}[\textit{i}]` must contain :math:`a_{\textit{i}}` (i.e., the coefficient of :math:`z^{{n-\textit{i}}}`), for :math:`\textit{i} = 0,1,\ldots,n`.
**scal** : bool, optional
Indicates whether or not the polynomial is to be scaled. See `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02agf.html#fcomments>`__ for advice on when it may be preferable to set :math:`\mathrm{scal} = \mathbf{False}` and for a description of the scaling strategy.
**Returns**
**z** : complex, ndarray, shape :math:`\left(n\right)`
The roots. Complex conjugate pairs of roots are stored in consecutive elements of :math:`\mathrm{z}`.
.. _c02ag-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{a}[0] = 0.0`.
Constraint: :math:`\mathrm{a}[0]\neq 0.0`.
(`errno` :math:`1`)
On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`.
Constraint: :math:`n \geq 1`.
(`errno` :math:`2`)
The iterative procedure has failed to converge. This error is very unlikely to occur. If it does, please contact `NAG <https://www.nag.com>`__ immediately, as some basic assumption for the arithmetic has been violated.
(`errno` :math:`3`)
``poly_real`` cannot evaluate :math:`P\left(z\right)` near some of its zeros without underflow. If this message occurs please contact `NAG <https://www.nag.com>`__.
(`errno` :math:`3`)
``poly_real`` cannot evaluate :math:`P\left(z\right)` near some of its zeros without overflow. If this message occurs please contact `NAG <https://www.nag.com>`__.
.. _c02ag-py2-py-notes:
**Notes**
`In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.`
``poly_real`` attempts to find all the roots of the :math:`n`\ th degree real polynomial equation
.. math::
P\left(z\right) = a_0z^n+a_1z^{{n-1}}+a_2z^{{n-2}} + \cdots +a_{{n-1}}z+a_n = 0\text{.}
The roots are located using a modified form of Laguerre's method, originally proposed by Smith (1967).
The method of Laguerre (see Wilkinson (1965)) can be described by the iterative scheme
.. math::
L\left(z_k\right) = z_{{k+1}}-z_k = \frac{{{-n}P\left(z_k\right)}}{{P^{\prime }\left(z_k\right)\pm \sqrt{H\left(z_k\right)}}}\text{,}
where :math:`H\left(z_k\right) = \left(n-1\right)\left[\left(n-1\right){\left(P^{\prime }\left(z_k\right)\right)}^2-nP\left(z_k\right)P^{{\prime \prime }}\left(z_k\right)\right]` and :math:`z_0` is specified.
The sign in the denominator is chosen so that the modulus of the Laguerre step at :math:`z_k`, viz. :math:`\left\lvert L\left(z_k\right)\right\rvert`, is as small as possible.
The method can be shown to be cubically convergent for isolated roots (real or complex) and linearly convergent for multiple roots.
The function generates a sequence of iterates :math:`z_1,z_2,z_3,\ldots \text{}`, such that :math:`\left\lvert P\left(z_{{k+1}}\right)\right\rvert < \left\lvert P\left(z_k\right)\right\rvert` and ensures that :math:`z_{{k+1}}+L\left(z_{{k+1}}\right)` 'roughly' lies inside a circular region of radius :math:`\left\lvert F\right\rvert` about :math:`z_k` known to contain a zero of :math:`P\left(z\right)`; that is, :math:`\left\lvert L\left(z_{{k+1}}\right)\right\rvert \leq \left\lvert F\right\rvert`, where :math:`F` denotes the Fejér bound (see Marden (1966)) at the point :math:`z_k`.
Following Smith (1967), :math:`F` is taken to be :math:`\mathrm{min}\left(B, {1.445nR}\right)`, where :math:`B` is an upper bound for the magnitude of the smallest zero given by
.. math::
B = 1.0001\times \mathrm{min}\left({\sqrt{n}L\left(z_k\right)}, \left\lvert r_1\right\rvert, \left\lvert a_n/a_0\right\rvert^{{1/n}}\right)\text{,}
:math:`r_1` is the zero :math:`X` of smaller magnitude of the quadratic equation
.. math::
\frac{{P^{{\prime \prime }}\left(z_k\right)}}{{2n\left(n-1\right)}}X^2+\frac{{P^{\prime }\left(z_k\right)}}{n}X+\frac{1}{2}P\left(z_k\right) = 0
and the Cauchy lower bound :math:`R` for the smallest zero is computed (using Newton's Method) as the positive root of the polynomial equation
.. math::
\left\lvert a_0\right\rvert z^n+\left\lvert a_1\right\rvert z^{{n-1}}+\left\lvert a_2\right\rvert z^{{n-2}} + \cdots +\left\lvert a_{{n-1}}\right\rvert z-\left\lvert a_n\right\rvert = 0\text{.}
Starting from the origin, successive iterates are generated according to the rule :math:`z_{{k+1}} = z_k+L\left(z_k\right)`, for :math:`k = 1,2,3,\ldots \text{}`, and :math:`L\left(z_k\right)` is 'adjusted' so that :math:`\left\lvert P\left(z_{{k+1}}\right)\right\rvert < \left\lvert P\left(z_k\right)\right\rvert` and :math:`\left\lvert L\left(z_{{k+1}}\right)\right\rvert \leq \left\lvert F\right\rvert`.
The iterative procedure terminates if :math:`P\left(z_{{k+1}}\right)` is smaller in absolute value than the bound on the rounding error in :math:`P\left(z_{{k+1}}\right)` and the current iterate :math:`z_p = z_{{k+1}}` is taken to be a zero of :math:`P\left(z\right)` (as is its conjugate :math:`\bar{z}_p` if :math:`z_p` is complex).
The deflated polynomial :math:`\tilde{P}\left(z\right) = P\left(z\right)/\left(z-z_p\right)` of degree :math:`n-1` if :math:`z_p` is real (:math:`\tilde{P}\left(z\right) = P\left(z\right)/\left(\left(z-z_p\right)\left(z-\bar{z}_p\right)\right)` of degree :math:`n-2` if :math:`z_p` is complex) is then formed, and the above procedure is repeated on the deflated polynomial until :math:`n < 3`, whereupon the remaining roots are obtained via the 'standard' closed formulae for a linear (:math:`n = 1`) or quadratic (:math:`n = 2`) equation.
Note that :meth:`quadratic_real`, :meth:`cubic_real` and :meth:`quartic_real` can be used to obtain the roots of a quadratic, cubic (:math:`n = 3`) and quartic (:math:`n = 4`) polynomial, respectively.
.. _c02ag-py2-py-references:
**References**
Marden, M, 1966, `Geometry of polynomials`, Mathematical Surveys (3), American Mathematical Society, Providence, RI
Smith, B T, 1967, `ZERPOL: a zero finding algorithm for polynomials using Laguerre's method`, Technical Report, Department of Computer Science, University of Toronto, Canada
Thompson, K W, 1991, `Error analysis for polynomial solvers`, Fortran Journal (Volume 3) (3), 10--13
Wilkinson, J H, 1965, `The Algebraic Eigenvalue Problem`, Oxford University Press, Oxford
"""
raise NotImplementedError
[docs]def quadratic_complex(a, b, c):
r"""
``quadratic_complex`` determines the roots of a quadratic equation with complex coefficients.
.. _c02ah-py2-py-doc:
For full information please refer to the NAG Library document for c02ah
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02ahf.html
.. _c02ah-py2-py-parameters:
**Parameters**
**a** : complex
:math:`a`, the coefficient of :math:`z^2`.
**b** : complex
:math:`b`, the coefficient of :math:`z`.
**c** : complex
:math:`c`, the constant coefficient.
**Returns**
**z_small** : complex
The smallest root in magnitude.
**z_large** : complex
The largest root in magnitude.
.. _c02ah-py2-py-errors:
**Warns**
**NagAlgorithmicWarning**
(`errno` :math:`1`)
On entry, :math:`\mathrm{a} = 0.0`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{a} = 0.0` and :math:`\mathrm{b} = 0.0`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{a} = 0.0` and the root :math:`{-\mathrm{c}}/\mathrm{b}` overflows: :math:`\mathrm{Re}\left(\mathrm{a}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Re}\left(\mathrm{c}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Re}\left(\mathrm{b}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Im}\left(\mathrm{a}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Im}\left(\mathrm{c}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Im}\left(\mathrm{b}\right) = \langle\mathit{\boldsymbol{value}}\rangle`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{c} = 0.0` and the root :math:`{-\mathrm{b}}/\mathrm{a}` overflows: :math:`\mathrm{Re}\left(\mathrm{c}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Re}\left(\mathrm{b}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Re}\left(\mathrm{a}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Im}\left(\mathrm{c}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Im}\left(\mathrm{b}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Im}\left(\mathrm{a}\right) = \langle\mathit{\boldsymbol{value}}\rangle`.
(`errno` :math:`5`)
On entry, :math:`B = \mathrm{max}\left(\left\lvert \mathrm{Re}\left(\mathrm{b}\right)\right\rvert, \left\lvert \mathrm{Im}\left(\mathrm{b}\right)\right\rvert \right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`A = \mathrm{max}\left(\left\lvert \mathrm{Re}\left(\mathrm{a}\right)\right\rvert, \left\lvert \mathrm{Im}\left(\mathrm{a}\right)\right\rvert \right) = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`C = \mathrm{max}\left(\left\lvert \mathrm{Re}\left(\mathrm{c}\right)\right\rvert, \left\lvert \mathrm{Im}\left(\mathrm{c}\right)\right\rvert \right) = \langle\mathit{\boldsymbol{value}}\rangle`. :math:`B` is so large that :math:`B^2` is indistinguishable from :math:`\left(B^2-4\times A\times C\right)` and the root :math:`{-\mathrm{b}}/\mathrm{a}` overflows. :math:`\mathrm{Re}\left(\mathrm{b}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Im}\left(\mathrm{b}\right) = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{Re}\left(\mathrm{a}\right) = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{Im}\left(\mathrm{a}\right) = \langle\mathit{\boldsymbol{value}}\rangle`.
.. _c02ah-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``quadratic_complex`` attempts to find the roots of the quadratic equation :math:`az^2+bz+c = 0` (where :math:`a`, :math:`b` and :math:`c` are complex coefficients), by carefully evaluating the 'standard' closed formula
.. math::
z = \frac{{-b\pm \sqrt{b^2-4ac}}}{{2a}}\text{.}
It is based on the function CQDRTC from Smith (1967).
**Note:** it is not necessary to scale the coefficients prior to calling the function.
.. _c02ah-py2-py-references:
**References**
Smith, B T, 1967, `ZERPOL: a zero finding algorithm for polynomials using Laguerre's method`, Technical Report, Department of Computer Science, University of Toronto, Canada
"""
raise NotImplementedError
[docs]def quadratic_real(a, b, c):
r"""
``quadratic_real`` determines the roots of a quadratic equation with real coefficients.
.. _c02aj-py2-py-doc:
For full information please refer to the NAG Library document for c02aj
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02ajf.html
.. _c02aj-py2-py-parameters:
**Parameters**
**a** : float
Must contain :math:`a`, the coefficient of :math:`z^2`.
**b** : float
Must contain :math:`b`, the coefficient of :math:`z`.
**c** : float
Must contain :math:`c`, the constant coefficient.
**Returns**
**z_small** : complex
The smallest root in magnitude.
**z_large** : complex
The largest root in magnitude.
.. _c02aj-py2-py-errors:
**Warns**
**NagAlgorithmicWarning**
(`errno` :math:`1`)
On entry, :math:`\mathrm{a} = 0.0`.
(`errno` :math:`2`)
On entry, :math:`\mathrm{a} = 0.0` and :math:`\mathrm{b} = 0.0`.
(`errno` :math:`3`)
On entry, :math:`\mathrm{a} = 0.0` and the root :math:`{-\mathrm{c}}/\mathrm{b}` overflows: :math:`\mathrm{a} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{c} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{b} = \langle\mathit{\boldsymbol{value}}\rangle`.
(`errno` :math:`4`)
On entry, :math:`\mathrm{c} = 0.0` and the root :math:`{-\mathrm{b}}/\mathrm{a}` overflows: :math:`\mathrm{c} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{b} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{a} = \langle\mathit{\boldsymbol{value}}\rangle`.
(`errno` :math:`5`)
On entry, :math:`\mathrm{b}` is so large that :math:`\mathrm{b}^2` is indistinguishable from :math:`\left(\mathrm{b}^2-4\times \mathrm{a}\times \mathrm{c}\right)` and the root :math:`{-\mathrm{b}}/\mathrm{a}` overflows: :math:`\mathrm{b} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{a} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{c} = \langle\mathit{\boldsymbol{value}}\rangle`.
.. _c02aj-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``quadratic_real`` attempts to find the roots of the quadratic equation :math:`az^2+bz+c = 0` (where :math:`a`, :math:`b` and :math:`c` are real coefficients), by carefully evaluating the 'standard' closed formula
.. math::
z = \frac{{-b\pm \sqrt{b^2-4ac}}}{{2a}}\text{.}
It is based on the function QDRTC from Smith (1967).
**Note:** it is not necessary to scale the coefficients prior to calling the function.
.. _c02aj-py2-py-references:
**References**
Smith, B T, 1967, `ZERPOL: a zero finding algorithm for polynomials using Laguerre's method`, Technical Report, Department of Computer Science, University of Toronto, Canada
"""
raise NotImplementedError
[docs]def cubic_real(u, r, s, t):
r"""
``cubic_real`` determines the roots of a cubic equation with real coefficients.
.. _c02ak-py2-py-doc:
For full information please refer to the NAG Library document for c02ak
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02akf.html
.. _c02ak-py2-py-parameters:
**Parameters**
**u** : float
:math:`u`, the coefficient of :math:`z^3`.
**r** : float
:math:`r`, the coefficient of :math:`z^2`.
**s** : float
:math:`s`, the coefficient of :math:`z`.
**t** : float
:math:`t`, the constant coefficient.
**Returns**
**zero** : complex, ndarray, shape :math:`\left(3\right)`
:math:`\mathrm{zero}[i-1]` contains the :math:`i`\ th root.
**errest** : float, ndarray, shape :math:`\left(3\right)`
:math:`\mathrm{errest}[i-1]` contains an approximate error estimate for the :math:`i`\ th root.
.. _c02ak-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{u} = 0.0`.
Constraint: :math:`\mathrm{u}\neq 0.0`.
(`errno` :math:`2`)
The companion matrix :math:`H` cannot be formed without overflow.
(`errno` :math:`3`)
Failure to converge in :meth:`lapackeig.dhseqr <naginterfaces.library.lapackeig.dhseqr>`.
.. _c02ak-py2-py-notes:
**Notes**
`In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.`
``cubic_real`` attempts to find the roots of the cubic equation
.. math::
uz^3+rz^2+sz+t = 0\text{,}
where :math:`u`, :math:`r`, :math:`s` and :math:`t` are real coefficients with :math:`u\neq 0`.
The roots are located by finding the eigenvalues of the associated :math:`3\times 3` (upper Hessenberg) companion matrix :math:`H` given by
.. math::
H = \begin{pmatrix}0&0&-t/u\\1&0&-s/u\\0&1&-r/u\end{pmatrix}\text{.}
The eigenvalues are obtained by a call to :meth:`lapackeig.dhseqr <naginterfaces.library.lapackeig.dhseqr>`.
Further details can be found in `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02akf.html#fcomments>`__.
To obtain the roots of a quartic equation, :meth:`quartic_real` can be used.
.. _c02ak-py2-py-references:
**References**
Golub, G H and Van Loan, C F, 1996, `Matrix Computations`, (3rd Edition), Johns Hopkins University Press, Baltimore
"""
raise NotImplementedError
[docs]def quartic_real(e, a, b, c, d):
r"""
``quartic_real`` determines the roots of a quartic equation with real coefficients.
.. _c02al-py2-py-doc:
For full information please refer to the NAG Library document for c02al
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02alf.html
.. _c02al-py2-py-parameters:
**Parameters**
**e** : float
:math:`e`, the coefficient of :math:`z^4`.
**a** : float
:math:`a`, the coefficient of :math:`z^3`.
**b** : float
:math:`b`, the coefficient of :math:`z^2`.
**c** : float
:math:`c`, the coefficient of :math:`z`.
**d** : float
:math:`d`, the constant coefficient.
**Returns**
**zero** : complex, ndarray, shape :math:`\left(4\right)`
:math:`\mathrm{zero}[i-1]` contains the :math:`i`\ th root.
**errest** : float, ndarray, shape :math:`\left(4\right)`
:math:`\mathrm{errest}[i-1]` contains an approximate error estimate for the :math:`i`\ th root.
.. _c02al-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{e} = 0.0`.
Constraint: :math:`\mathrm{e}\neq 0.0`.
(`errno` :math:`2`)
The companion matrix :math:`H` cannot be formed without overflow.
(`errno` :math:`3`)
Failure to converge in :meth:`lapackeig.dhseqr <naginterfaces.library.lapackeig.dhseqr>`.
.. _c02al-py2-py-notes:
**Notes**
`In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.`
``quartic_real`` attempts to find the roots of the quartic equation
.. math::
ez^4+az^3+bz^2+cz+d = 0\text{,}
where :math:`e`, :math:`a`, :math:`b`, :math:`c` and :math:`d` are real coefficients with :math:`e\neq 0`.
The roots are located by finding the eigenvalues of the associated :math:`4\times 4` (upper Hessenberg) companion matrix :math:`H` given by
.. math::
H = \begin{pmatrix}0&0&0&-d/e\\1&0&0&-c/e\\0&1&0&-b/e\\0&0&1&-a/e\end{pmatrix}\text{.}
The eigenvalues are obtained by a call to :meth:`lapackeig.dhseqr <naginterfaces.library.lapackeig.dhseqr>`.
Further details can be found in `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02alf.html#fcomments>`__.
To obtain the roots of a cubic equation, :meth:`cubic_real` can be used.
.. _c02al-py2-py-references:
**References**
Golub, G H and Van Loan, C F, 1996, `Matrix Computations`, (3rd Edition), Johns Hopkins University Press, Baltimore
"""
raise NotImplementedError
[docs]def cubic_complex(u, r, s, t):
r"""
``cubic_complex`` determines the roots of a cubic equation with complex coefficients.
.. _c02am-py2-py-doc:
For full information please refer to the NAG Library document for c02am
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02amf.html
.. _c02am-py2-py-parameters:
**Parameters**
**u** : complex
:math:`u`, the coefficient of :math:`z^3`.
**r** : complex
:math:`r`, the coefficient of :math:`z^2`.
**s** : complex
:math:`s`, the coefficient of :math:`z`.
**t** : complex
:math:`t`, the constant coefficient.
**Returns**
**zero** : complex, ndarray, shape :math:`\left(3\right)`
:math:`\mathrm{zero}[i-1]` contains the :math:`i`\ th root.
**errest** : float, ndarray, shape :math:`\left(3\right)`
:math:`\mathrm{errest}[i-1]` contains an approximate error estimate for the :math:`i`\ th root.
.. _c02am-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{u} = \left(0.0, 0.0\right)`.
Constraint: :math:`\mathrm{u}\neq \left(0.0, 0.0\right)`
(`errno` :math:`2`)
The companion matrix :math:`H` cannot be formed without overflow.
(`errno` :math:`3`)
Failure to converge in :meth:`lapackeig.zhseqr <naginterfaces.library.lapackeig.zhseqr>`.
.. _c02am-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``cubic_complex`` attempts to find the roots of the cubic equation
.. math::
uz^3+rz^2+sz+t = 0\text{,}
where :math:`u`, :math:`r`, :math:`s` and :math:`t` are complex coefficients with :math:`u\neq 0`.
The roots are located by finding the eigenvalues of the associated :math:`3\times 3` (upper Hessenberg) companion matrix :math:`H` given by
.. math::
H = \begin{pmatrix}0&0&-t/u\\1&0&-s/u\\0&1&-r/u\end{pmatrix}\text{.}
The eigenvalues are obtained by a call to :meth:`lapackeig.zhseqr <naginterfaces.library.lapackeig.zhseqr>`.
Further details can be found in `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02amf.html#fcomments>`__.
To obtain the roots of a quadratic equation, :meth:`quadratic_complex` can be used.
.. _c02am-py2-py-references:
**References**
Golub, G H and Van Loan, C F, 1996, `Matrix Computations`, (3rd Edition), Johns Hopkins University Press, Baltimore
"""
raise NotImplementedError
[docs]def quartic_complex(e, a, b, c, d):
r"""
``quartic_complex`` determines the roots of a quartic equation with complex coefficients.
.. _c02an-py2-py-doc:
For full information please refer to the NAG Library document for c02an
https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02anf.html
.. _c02an-py2-py-parameters:
**Parameters**
**e** : complex
:math:`e`, the coefficient of :math:`z^4`.
**a** : complex
:math:`a`, the coefficient of :math:`z^3`.
**b** : complex
:math:`b`, the coefficient of :math:`z^2`.
**c** : complex
:math:`c`, the coefficient of :math:`z`.
**d** : complex
:math:`d`, the constant coefficient.
**Returns**
**zero** : complex, ndarray, shape :math:`\left(4\right)`
:math:`\mathrm{zero}[i-1]` contains the :math:`i`\ th root.
**errest** : float, ndarray, shape :math:`\left(4\right)`
:math:`\mathrm{errest}[i-1]` contains an approximate error estimate for the :math:`i`\ th root.
.. _c02an-py2-py-errors:
**Raises**
**NagValueError**
(`errno` :math:`1`)
On entry, :math:`\mathrm{e} = \left(0.0, 0.0\right)`.
Constraint: :math:`\mathrm{e}\neq \left(0.0, 0.0\right)`
(`errno` :math:`2`)
The companion matrix :math:`H` cannot be formed without overflow.
(`errno` :math:`3`)
Failure to converge in :meth:`lapackeig.zhseqr <naginterfaces.library.lapackeig.zhseqr>`.
.. _c02an-py2-py-notes:
**Notes**
`No equivalent traditional C interface for this routine exists in the NAG Library.`
``quartic_complex`` attempts to find the roots of the quartic equation
.. math::
ez^4+az^3+bz^2+cz+d = 0\text{,}
where :math:`e`, :math:`a`, :math:`b`, :math:`c` and :math:`d` are complex coefficients with :math:`e\neq 0`.
The roots are located by finding the eigenvalues of the associated :math:`4\times 4` (upper Hessenberg) companion matrix :math:`H` given by
.. math::
H = \begin{pmatrix}0&0&0&-d/e\\1&0&0&-c/e\\0&1&0&-b/e\\0&0&1&-a/e\end{pmatrix}\text{.}
The eigenvalues are obtained by a call to :meth:`lapackeig.zhseqr <naginterfaces.library.lapackeig.zhseqr>`.
Further details can be found in `Further Comments <https://support.nag.com/numeric/nl/nagdoc_30/flhtml/c02/c02anf.html#fcomments>`__.
To obtain the roots of a cubic equation, :meth:`cubic_complex` can be used.
.. _c02an-py2-py-references:
**References**
Golub, G H and Van Loan, C F, 1996, `Matrix Computations`, (3rd Edition), Johns Hopkins University Press, Baltimore
"""
raise NotImplementedError
