Polynomials revision c981c48f5bc9aefeffc0bcb0cc3934c2fae179dd 
1#ifndef EIGEN_POLYNOMIALS_MODULE_H
2#define EIGEN_POLYNOMIALS_MODULE_H
3
4#include <Eigen/Core>
5
6#include <Eigen/src/Core/util/DisableStupidWarnings.h>
7
8#include <Eigen/Eigenvalues>
9
10// Note that EIGEN_HIDE_HEAVY_CODE has to be defined per module
11#if (defined EIGEN_EXTERN_INSTANTIATIONS) && (EIGEN_EXTERN_INSTANTIATIONS>=2)
12  #ifndef EIGEN_HIDE_HEAVY_CODE
13  #define EIGEN_HIDE_HEAVY_CODE
14  #endif
15#elif defined EIGEN_HIDE_HEAVY_CODE
16  #undef EIGEN_HIDE_HEAVY_CODE
17#endif
18
19/** \ingroup Unsupported_modules
20  * \defgroup Polynomials_Module Polynomials module
21  *
22  *
23  *
24  * \brief This module provides a QR based polynomial solver.
25	*
26  * To use this module, add
27  * \code
28  * #include <unsupported/Eigen/Polynomials>
29  * \endcode
30	* at the start of your source file.
31  */
32
33#include "src/Polynomials/PolynomialUtils.h"
34#include "src/Polynomials/Companion.h"
35#include "src/Polynomials/PolynomialSolver.h"
36
37/**
38	\page polynomials Polynomials defines functions for dealing with polynomials
39	and a QR based polynomial solver.
40	\ingroup Polynomials_Module
41
42	The remainder of the page documents first the functions for evaluating, computing
43	polynomials, computing estimates about polynomials and next the QR based polynomial
44	solver.
45
46	\section polynomialUtils convenient functions to deal with polynomials
47	\subsection roots_to_monicPolynomial
48	The function
49	\code
50	void roots_to_monicPolynomial( const RootVector& rv, Polynomial& poly )
51	\endcode
52	computes the coefficients \f$ a_i \f$ of
53
54	\f$ p(x) = a_0 + a_{1}x + ... + a_{n-1}x^{n-1} + x^n \f$
55
56	where \f$ p \f$ is known through its roots i.e. \f$ p(x) = (x-r_1)(x-r_2)...(x-r_n) \f$.
57
58	\subsection poly_eval
59	The function
60	\code
61	T poly_eval( const Polynomials& poly, const T& x )
62	\endcode
63	evaluates a polynomial at a given point using stabilized H&ouml;rner method.
64
65	The following code: first computes the coefficients in the monomial basis of the monic polynomial that has the provided roots;
66	then, it evaluates the computed polynomial, using a stabilized H&ouml;rner method.
67
68	\include PolynomialUtils1.cpp
69  Output: \verbinclude PolynomialUtils1.out
70
71	\subsection Cauchy bounds
72	The function
73	\code
74	Real cauchy_max_bound( const Polynomial& poly )
75	\endcode
76	provides a maximum bound (the Cauchy one: \f$C(p)\f$) for the absolute value of a root of the given polynomial i.e.
77	\f$ \forall r_i \f$ root of \f$ p(x) = \sum_{k=0}^d a_k x^k \f$,
78	\f$ |r_i| \le C(p) = \sum_{k=0}^{d} \left | \frac{a_k}{a_d} \right | \f$
79	The leading coefficient \f$ p \f$: should be non zero \f$a_d \neq 0\f$.
80
81
82	The function
83	\code
84	Real cauchy_min_bound( const Polynomial& poly )
85	\endcode
86	provides a minimum bound (the Cauchy one: \f$c(p)\f$) for the absolute value of a non zero root of the given polynomial i.e.
87	\f$ \forall r_i \neq 0 \f$ root of \f$ p(x) = \sum_{k=0}^d a_k x^k \f$,
88	\f$ |r_i| \ge c(p) = \left( \sum_{k=0}^{d} \left | \frac{a_k}{a_0} \right | \right)^{-1} \f$
89
90
91
92
93	\section QR polynomial solver class
94	Computes the complex roots of a polynomial by computing the eigenvalues of the associated companion matrix with the QR algorithm.
95	
96	The roots of \f$ p(x) = a_0 + a_1 x + a_2 x^2 + a_{3} x^3 + x^4 \f$ are the eigenvalues of
97	\f$
98	\left [
99	\begin{array}{cccc}
100	0 & 0 &  0 & a_0 \\
101	1 & 0 &  0 & a_1 \\
102	0 & 1 &  0 & a_2 \\
103	0 & 0 &  1 & a_3
104	\end{array} \right ]
105	\f$
106
107	However, the QR algorithm is not guaranteed to converge when there are several eigenvalues with same modulus.
108
109	Therefore the current polynomial solver is guaranteed to provide a correct result only when the complex roots \f$r_1,r_2,...,r_d\f$ have distinct moduli i.e.
110	
111	\f$ \forall i,j \in [1;d],~ \| r_i \| \neq \| r_j \| \f$.
112
113	With 32bit (float) floating types this problem shows up frequently.
114  However, almost always, correct accuracy is reached even in these cases for 64bit
115  (double) floating types and small polynomial degree (<20).
116
117	\include PolynomialSolver1.cpp
118	
119	In the above example:
120	
121	-# a simple use of the polynomial solver is shown;
122	-# the accuracy problem with the QR algorithm is presented: a polynomial with almost conjugate roots is provided to the solver.
123	Those roots have almost same module therefore the QR algorithm failed to converge: the accuracy
124	of the last root is bad;
125	-# a simple way to circumvent the problem is shown: use doubles instead of floats.
126
127  Output: \verbinclude PolynomialSolver1.out
128*/
129
130#include <Eigen/src/Core/util/ReenableStupidWarnings.h>
131
132#endif // EIGEN_POLYNOMIALS_MODULE_H
133/* vim: set filetype=cpp et sw=2 ts=2 ai: */
134