6#ifndef THEORETICA_POLYNOMIAL_H
7#define THEORETICA_POLYNOMIAL_H
9#ifndef THEORETICA_NO_PRINT
14#include "../core/real_analysis.h"
15#include "../algebra/vec.h"
16#include "../complex/complex.h"
17#include "../complex/complex_analysis.h"
24 template<
typename Type = real>
27 std::vector<Type> coeff;
70 template<
typename EvalType = Type>
79 for (
unsigned int i = 0; i < coeff.size(); ++i)
80 sum = coeff[coeff.size() - i - 1] + x *
sum;
87 template<
typename EvalType = Type>
97 for (
int i = coeff.size() - 1; i >= 0; --i) {
109 for (
unsigned int i = coeff.size() - 1; i >= 0; --i) {
128 r.coeff.resize(
max(coeff.size(),
p.size()));
132 for (; i <
min(coeff.size(),
p.size()); ++i)
133 r[i] = coeff[i] +
p[i];
135 if (coeff.size() >
p.size())
136 for (; i < coeff.size(); ++i)
139 for (; i <
p.size(); ++i)
150 r.coeff.resize(
max(coeff.size(),
p.size()));
154 for (; i <
min(
r.size(),
p.size()); ++i)
155 r[i] = coeff[i] -
p.get(i);
157 if (coeff.size() >
p.size())
158 for (; i < coeff.size(); ++i)
161 for (; i <
p.size(); ++i)
172 r.coeff.resize(this->
size() + p.
size() - 1);
174 for (
unsigned int i = 0; i <
size(); ++i) {
175 for (
unsigned int j = 0;
j <
p.size(); ++
j) {
176 r[i +
j] += coeff[i] *
p.get(
j);
205 const unsigned int r_order =
r.find_order();
228 TH_MATH_ERROR(
"polynomial::operator/", i, NO_ALGO_CONVERGENCE);
241 for (
unsigned int i = 0; i <
size(); ++i)
257 for (
unsigned int i = 0; i <
size(); ++i)
268 if(coeff.size() <
p.size())
269 coeff.resize(
p.size(),
Type(0));
271 for (
unsigned int i = 0; i <
min(
size(),
p.size()); ++i)
272 coeff[i] +=
p.
get(i);
282 if(coeff.size() <
p.size())
283 coeff.resize(
p.size(),
Type(0));
285 for (
unsigned int i = 0; i <
min(coeff.size(),
p.size()); ++i)
286 coeff[i] -=
p.
get(i);
298 for (
unsigned int i = 0; i < coeff.size(); ++i)
299 for (
unsigned int j = 0;
j <
p.size(); ++
j)
300 r[i +
j] += coeff[i] *
p.get(
j);
310 for (
unsigned int i = 0; i < coeff.size(); ++i)
319 return (*
this = (*
this / a));
331 for (
unsigned int i = 0; i < coeff.size(); ++i)
342 for (
unsigned int i = 0; i <
n; ++i) {
348 for (
unsigned int i = 0; i <
size(); ++i) {
353 for (
unsigned int i = 0; i <
other.size(); ++i) {
366 return coeff.begin();
380 return coeff.cbegin();
402 const Type p = coeff[1] / coeff[2];
403 const Type q = coeff[0] / coeff[2];
431 const std::vector<Type>&
roots) {
436 for (
unsigned int i = 0; i <
roots.size(); ++i)
447 m.coeff = std::vector<Type>(
order + 1,
Type(0));
475#ifndef THEORETICA_NO_PRINT
479 const std::string&
unknown =
"x",
482 std::stringstream
res;
484 const int sz = coeff.size();
486 for (
int i =
sz - 1; i >= 0; --i) {
491 res << (coeff[i] >= 0 ?
"+ " :
"- ");
510 inline operator std::string() {
518 return out <<
obj.to_string();
A polynomial of arbitrary order.
Definition polynomial.h:25
polynomial & operator/=(const polynomial &a)
Multiply a polynomial by a scalar value.
Definition polynomial.h:318
polynomial & operator+=(const polynomial &p)
Sum a polynomial to this one.
Definition polynomial.h:265
~polynomial()
Default destructor.
Definition polynomial.h:42
auto begin()
Get an iterator for the first coefficient of the polynomial.
Definition polynomial.h:365
polynomial(std::initializer_list< Type > l)
Initialize from an std::initializer_list.
Definition polynomial.h:39
polynomial(const std::vector< Type > &c)
Initialize from an std::vector.
Definition polynomial.h:36
auto cbegin()
Get a const iterator for the first coefficient of the polynomial.
Definition polynomial.h:379
EvalType eval(EvalType x) const
Evaluate the polynomial using x as variable.
Definition polynomial.h:71
polynomial operator/(const polynomial &d) const
Polynomial division.
Definition polynomial.h:185
void trim()
Remove trailing zero coefficients.
Definition polynomial.h:107
Type get(int i) const
Get the i-th by value.
Definition polynomial.h:52
EvalType operator()(EvalType x) const
Evaluate the polynomial using x as variable.
Definition polynomial.h:88
polynomial & operator*=(const polynomial &p)
Multiply two polynomials.
Definition polynomial.h:293
polynomial(Type a)
Initialize as a constant.
Definition polynomial.h:33
auto end()
Get an iterator for the one plus last coefficient of the polynomial.
Definition polynomial.h:372
vec< complex<>, 2 > quadratic_roots() const
Compute the roots of a quadratic polynomial.
Definition polynomial.h:392
friend std::ostream & operator<<(std::ostream &out, const polynomial &obj)
Stream the polynomial in string representation to an output stream (std::ostream)
Definition polynomial.h:517
polynomial & operator-=(const polynomial &p)
Subtract a polynomial from this one.
Definition polynomial.h:279
const Type & operator[](unsigned int i) const
Return the nth order coefficient.
Definition polynomial.h:58
static polynomial< Type > monomial(Type c, unsigned int order)
Returns a monomial of the given degree and coefficient.
Definition polynomial.h:444
bool operator==(const polynomial &other) const
Check whether two polynomials are equal.
Definition polynomial.h:339
polynomial operator*(const polynomial &p) const
Multiply two polynomials.
Definition polynomial.h:169
polynomial & operator/=(Type a)
Divide a polynomial by a scalar value.
Definition polynomial.h:324
Type & operator[](unsigned int i)
Return the nth order coefficient.
Definition polynomial.h:64
polynomial()
Initialize as an empty polynomial.
Definition polynomial.h:30
unsigned int find_order() const
Find the true order of the polynomial (ignoring trailing null coefficients)
Definition polynomial.h:95
auto cend()
Get a const iterator for the one plus last coefficient of the polynomial.
Definition polynomial.h:386
size_t size() const
Get the number of coefficients.
Definition polynomial.h:119
polynomial operator+(const polynomial &p) const
Sum two polynomials.
Definition polynomial.h:125
polynomial operator-(const polynomial &p) const
Subtract a polynomial from another.
Definition polynomial.h:147
polynomial operator/(Type a) const
Divide a polynomial by a scalar.
Definition polynomial.h:249
std::string to_string(const std::string &unknown="x", const std::string &exponentiation="^") const
Convert the polynomial to string representation.
Definition polynomial.h:478
polynomial & operator*=(Type a)
Multiply a polynomial by a scalar value.
Definition polynomial.h:308
polynomial operator*(Type a) const
Multiply a polynomial by a scalar.
Definition polynomial.h:237
static polynomial< Type > from_roots(const std::vector< Type > &roots)
Construct a polynomial from its roots.
Definition polynomial.h:430
Type & at(int i)
Access i-th coefficient.
Definition polynomial.h:46
#define TH_MATH_ERROR(F_NAME, VALUE, EXCEPTION)
TH_MATH_ERROR is a macro which throws exceptions or modifies errno (depending on which compiling opti...
Definition error.h:219
Main namespace of the library which contains all functions and objects.
Definition algebra.h:27
auto min(const Vector &X)
Finds the minimum value inside a dataset.
Definition dataset.h:351
dual2 sqrt(dual2 x)
Compute the square root of a second order dual number.
Definition dual2_functions.h:54
dual2 abs(dual2 x)
Compute the absolute value of a second order dual number.
Definition dual2_functions.h:198
std::remove_reference_t< decltype(std::declval< Structure >()[0])> vector_element_t
Extract the type of a vector (or any indexable container) from its operator[].
Definition core_traits.h:134
auto sum(const Vector &X)
Compute the sum of a vector of real values using pairwise summation to reduce round-off error.
Definition dataset.h:219
std::vector< real > roots(RealFunction f, real a, real b, real tol=OPTIMIZATION_TOL, real steps=10)
Find the roots of a univariate real function inside a given interval, by first searching for candidat...
Definition roots.h:650
auto max(const Vector &X)
Finds the maximum value inside a dataset.
Definition dataset.h:330
constexpr real MACH_EPSILON
Machine epsilon for the real type.
Definition constants.h:207
int sgn(real x)
Return the sign of x (1 if positive, -1 if negative, 0 if null)
Definition real_analysis.h:259
dual2 square(dual2 x)
Return the square of a second order dual number.
Definition dual2_functions.h:23
real nan()
Return a quiet NaN number in floating point representation.
Definition error.h:54