1 #ifndef DOUBLE_POLY_H_ 2 #define DOUBLE_POLY_H_ 3 4 // IWYU pragma: no_include "mpz_poly.h" 5 // the fwd-decl is enough 6 7 #include <stdint.h> // for uint32_t 8 #include <stdio.h> // for FILE 9 #include "macros.h" // for GNUC_VERSION_ATLEAST 10 #ifdef __cplusplus 11 #include <string> // for string 12 #endif 13 14 #ifdef DOUBLE_POLY_EXPOSE_COMPLEX_FUNCTIONS 15 #include <complex.h> // IWYU pragma: keep 16 #endif 17 18 #ifndef MPZ_POLY_H_ 19 typedef struct mpz_poly_s * mpz_poly_ptr; 20 typedef const struct mpz_poly_s * mpz_poly_srcptr; 21 #endif 22 23 /* floating point polynomials */ 24 25 #ifdef __cplusplus 26 extern "C" { 27 #endif 28 29 struct double_poly_s { 30 int alloc; 31 int deg; 32 double *coeff; 33 }; 34 35 typedef struct double_poly_s double_poly[1]; 36 #ifndef MPZ_POLY_H_ 37 /* mpz_poly.h forward-declares these. Don't do it twice */ 38 typedef struct double_poly_s * double_poly_ptr; 39 typedef const struct double_poly_s * double_poly_srcptr; 40 #endif 41 42 /* double_poly.c */ 43 void double_poly_init (double_poly_ptr, int deg); 44 void double_poly_realloc (double_poly_ptr, int nc); 45 void double_poly_clear (double_poly_ptr); 46 void double_poly_swap(double_poly_ptr p, double_poly_ptr q); 47 void double_poly_set (double_poly_ptr, double_poly_srcptr); 48 49 void double_poly_set_zero (double_poly_ptr r); 50 void double_poly_set_xi(double_poly_ptr s, int i); 51 52 void double_poly_cleandeg(double_poly_ptr f, int deg); 53 54 int double_poly_cmp(double_poly_srcptr a, double_poly_srcptr b); 55 56 double double_poly_eval (double_poly_srcptr, double); 57 double double_poly_eval_homogeneous (double_poly_srcptr p, double x, double y); 58 double double_poly_eval_safe (double_poly_srcptr, double); 59 double double_poly_dichotomy (double_poly_srcptr, double, double, double, 60 unsigned int); 61 void double_poly_derivative(double_poly_ptr, double_poly_srcptr); 62 void double_poly_mul(double_poly_ptr, double_poly_srcptr, double_poly_srcptr); 63 void double_poly_add(double_poly_ptr, double_poly_srcptr, double_poly_srcptr); 64 void double_poly_sub(double_poly_ptr, double_poly_srcptr, double_poly_srcptr); 65 void double_poly_addmul(double_poly_ptr, double_poly_srcptr, double_poly_srcptr); 66 void double_poly_submul(double_poly_ptr, double_poly_srcptr, double_poly_srcptr); 67 void double_poly_addmul_double(double_poly_ptr, double_poly_srcptr, double); 68 void double_poly_submul_double(double_poly_ptr, double_poly_srcptr, double); 69 void double_poly_neg(double_poly_ptr, double_poly_srcptr); 70 void double_poly_revert(double_poly_ptr, double_poly_srcptr); 71 double double_poly_bound_roots (double_poly_srcptr p); 72 unsigned int double_poly_compute_roots(double *, double_poly_srcptr, double); 73 unsigned int double_poly_compute_all_roots_with_bound(double *, 74 double_poly_srcptr, 75 double); 76 unsigned int double_poly_compute_all_roots(double *, double_poly_srcptr); 77 void double_poly_print (FILE *, double_poly_srcptr, char *variable_name); 78 int double_poly_asprint (char **t, double_poly_srcptr p, char *variable_name); 79 void double_poly_set_mpz_poly (double_poly_ptr p, mpz_poly_srcptr q); 80 81 double double_poly_resultant(double_poly_srcptr p, double_poly_srcptr q); 82 void double_poly_mul_double(double_poly_ptr f, double_poly_srcptr g, 83 double mul); 84 double double_poly_div_linear(double_poly_ptr q, double_poly_srcptr p, const double r); 85 void double_poly_set_string(double_poly_ptr poly, const char *str); 86 87 #ifdef DOUBLE_POLY_EXPOSE_COMPLEX_FUNCTIONS 88 /* we use _Complex and not complex below, because it's mandatory in order 89 * to be nice to C++ code */ 90 void double_poly_complex_roots(double _Complex *roots, double_poly_srcptr); 91 void double_poly_complex_roots_long(long double _Complex *roots, double_poly_srcptr); 92 uint32_t poly_roots_double(double *poly, uint32_t degree, double _Complex *roots); 93 uint32_t poly_roots_longdouble(double *poly, uint32_t degree, long double _Complex *roots); 94 #endif 95 96 #ifdef __cplusplus 97 } 98 #endif 99 100 #ifdef __cplusplus 101 /* This is sort of a generic way to write a c++ equivalent to the C type. 102 * The first-class citizen in the cado-nfs code is (still) the C type, so 103 * we're definitely bound to have a few infelicities here: 104 * - the type name can't be the same because of the size-1 array trick 105 * in C. 106 * - the C type is embedded as a member x for the same reason. 107 * - most operations on the C type should go through the member x 108 * (however, the conversions we have to _ptr and _srcptr can ease 109 * things a bit). 110 */ 111 struct cxx_double_poly { 112 double_poly x; 113 cxx_double_poly(int deg = -1) { double_poly_init(x, deg); } cxx_double_polycxx_double_poly114 cxx_double_poly(double_poly_srcptr f) { double_poly_init(x, -1); double_poly_set(x, f); } ~cxx_double_polycxx_double_poly115 ~cxx_double_poly() { double_poly_clear(x); } cxx_double_polycxx_double_poly116 cxx_double_poly(cxx_double_poly const & o) { 117 double_poly_init(x, -1); 118 double_poly_set(x, o.x); 119 } 120 cxx_double_poly & operator=(cxx_double_poly const & o) { 121 double_poly_set(x, o.x); 122 return *this; 123 } 124 #if __cplusplus >= 201103L cxx_double_polycxx_double_poly125 cxx_double_poly(cxx_double_poly && o) { 126 double_poly_init(x, -1); 127 double_poly_swap(x, o.x); 128 } 129 cxx_double_poly& operator=(cxx_double_poly && o) { 130 double_poly_swap(x, o.x); 131 return *this; 132 } 133 #endif double_poly_ptrcxx_double_poly134 operator double_poly_ptr() { return x; } double_poly_srcptrcxx_double_poly135 operator double_poly_srcptr() const { return x; } 136 double_poly_ptr operator->() { return x; } 137 double_poly_srcptr operator->() const { return x; } 138 std::string print_poly(std::string const& var) const; 139 }; 140 #if GNUC_VERSION_ATLEAST(4,3,0) 141 extern void double_poly_init(cxx_double_poly & pl, int) __attribute__((error("double_poly_init must not be called on a double_poly reference -- it is the caller's business (via a ctor)"))); 142 extern void double_poly_clear(cxx_double_poly & pl) __attribute__((error("double_poly_clear must not be called on a double_poly reference -- it is the caller's business (via a dtor)"))); 143 #endif 144 145 #endif 146 147 148 #endif /* DOUBLE_POLY_H_ */ 149