FastPolyEval  1.0
Fast Evaluation of Real and Complex Polynomials
poly.h
Go to the documentation of this file.
1 //
2 // poly.h
3 //
4 // Authors: Nicolae Mihalache & François Vigneron
5 //
6 // This software is released under BSD licence, with an attribution clause (see Licence file).
7 // Please cite the reference below if you use or distribute this software.
8 //
9 // • [1] R. Anton, N. Mihalache & F. Vigneron. Fast evaluation of real and complex polynomials. 2022.
10 // https://hal.archives-ouvertes.fr/hal-03820369
11 //
12 // Copyright 2022 Univ. Paris-Est Créteil & Univ. de Reims Champagne-Ardenne.//
13 
19 #ifndef poly_h
20 #define poly_h
21 
22 #include <mpfr.h>
23 
24 #include "ntypes.h"
25 #include "comp.h"
26 
27 // MARK: data types
28 
31 typedef struct {
35  bool modified;
36  mpfr_t buf1;
37  mpfr_t buf2;
38 } poly_struct;
39 
43 typedef poly_struct poly_t[1];
44 
46 typedef poly_struct *poly;
47 
48 // MARK: functions prototypes
49 
58 poly poly_new(deg_t degree, prec_t prec);
59 
69 poly poly_from_roots(comp_ptr roots, deg_t degree, prec_t prec);
70 
77 bool poly_free(poly P);
78 
86 bool poly_set(poly P, comp coeff, deg_t ind);
87 
95 bool poly_eval(comp res, poly P, comp z);
96 
104 bool poly_eval_r(comp res, poly P, mpfr_t x);
105 
112 
120 
128 
136 
143 
144 #endif /* poly_h */
Definition of MPFR complex numbers.
comp_struct comp[1]
Practical wrapper for comp_struct.
Definition: comp.h:39
Definition of basic types.
ulong deg_t
The integer number type to use for polynomial degrees and indexes.
Definition: ntypes.h:128
ulong prec_t
The integer number type to use for polynomial degrees and indexes.
Definition: ntypes.h:188
poly poly_derivative(poly P)
Computes the derivative of P.
poly poly_sqr(poly P)
Computes the square of P.
bool poly_eval(comp res, poly P, comp z)
Evaluates P(z) using Horner's method.
bool poly_set(poly P, comp coeff, deg_t ind)
Sets the coefficient of the polynomial P corresponding to the power ind to coeff.
poly poly_diff(poly P, poly Q)
Computes P-Q.
poly_struct * poly
Convenience pointer to poly_struct.
Definition: poly.h:46
bool poly_eval_r(comp res, poly P, mpfr_t x)
Evaluates P(x) using Horner's method.
poly poly_from_roots(comp_ptr roots, deg_t degree, prec_t prec)
Returns a new complex polynomial given the list of its roots, with coefficients of precision prec.
poly poly_prod(poly P, poly Q)
Computes P*Q.
poly poly_sum(poly P, poly Q)
Computes P+Q.
poly_struct poly_t[1]
Practical wrapper for poly_struct.
Definition: poly.h:43
poly poly_new(deg_t degree, prec_t prec)
Returns a new complex polynomial of given degree, with coefficients of precision prec.
bool poly_free(poly P)
Frees all the memory used by the polynomial P, assuming the struct has been allocated with malloc(),...
Multi-precision floating point complex numbers.
Definition: comp.h:31
Polynomial with multi-precision floating point complex coefficients.
Definition: poly.h:31
prec_t prec
the precision of the coefficients, in bits
Definition: poly.h:33
bool modified
the status of the coefficients
Definition: poly.h:35
mpfr_t buf2
another buffer
Definition: poly.h:37
comp_ptr a
the coefficients
Definition: poly.h:34
mpfr_t buf1
a buffer
Definition: poly.h:36
deg_t degree
the degree of the polynomial
Definition: poly.h:32