1 #pragma once
2 /* #[ License : */
3 /*
4 * Copyright (C) 1984-2017 J.A.M. Vermaseren
5 * When using this file you are requested to refer to the publication
6 * J.A.M.Vermaseren "New features of FORM" math-ph/0010025
7 * This is considered a matter of courtesy as the development was paid
8 * for by FOM the Dutch physics granting agency and we would like to
9 * be able to track its scientific use to convince FOM of its value
10 * for the community.
11 *
12 * This file is part of FORM.
13 *
14 * FORM is free software: you can redistribute it and/or modify it under the
15 * terms of the GNU General Public License as published by the Free Software
16 * Foundation, either version 3 of the License, or (at your option) any later
17 * version.
18 *
19 * FORM is distributed in the hope that it will be useful, but WITHOUT ANY
20 * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
21 * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
22 * details.
23 *
24 * You should have received a copy of the GNU General Public License along
25 * with FORM. If not, see <http://www.gnu.org/licenses/>.
26 */
27 /* #] License : */
28
29 extern "C" {
30 #include "form3.h"
31 }
32
33 #include <string>
34 #include <vector>
35
36 // macros for tform
37 #ifndef WITHPTHREADS
38 #define POLY_GETIDENTITY(X)
39 #define POLY_STOREIDENTITY
40 #else
41 #define POLY_GETIDENTITY(X) ALLPRIVATES *B = (X).Bpointer
42 #define POLY_STOREIDENTITY Bpointer = B
43 #endif
44
45 // maximum size of the hash table used for multiplication and division
46
47 const int POLY_MAX_HASH_SIZE = MiN(1<<20, MAXPOSITIVE);
48
49 class poly {
50
51 public:
52
53 // variables
54 #ifdef WITHPTHREADS
55 ALLPRIVATES *Bpointer;
56 #endif
57
58 WORD *terms;
59 LONG size_of_terms;
60 WORD modp,modn;
61
62 // constructors/destructor
63 poly (PHEAD int, WORD=-1, WORD=1);
64 poly (PHEAD const UWORD *, WORD, WORD=-1, WORD=1);
65 poly (const poly &, WORD=-1, WORD=1);
66 ~poly ();
67
68 // operators
69 poly& operator+= (const poly&);
70 poly& operator-= (const poly&);
71 poly& operator*= (const poly&);
72 poly& operator/= (const poly&);
73 poly& operator%= (const poly&);
74
75 const poly operator+ (const poly&) const;
76 const poly operator- (const poly&) const;
77 const poly operator* (const poly&) const;
78 const poly operator/ (const poly&) const;
79 const poly operator% (const poly&) const;
80
81 bool operator== (const poly&) const;
82 bool operator!= (const poly&) const;
83 poly& operator= (const poly &);
84 WORD& operator[] (int);
85 const WORD& operator[] (int) const;
86
87 // memory management
88 void termscopy (const WORD *, int, int);
89 void check_memory(int);
90 void expand_memory(int);
91
92 // check type of polynomial
93 bool is_zero () const;
94 bool is_one () const;
95 bool is_integer () const;
96 bool is_monomial () const;
97 int is_dense_univariate () const;
98
99 // properties
100 int sign () const;
101 int degree (int) const;
102 int total_degree () const;
103 int first_variable () const;
104 int number_of_terms () const;
105 const std::vector<int> all_variables () const;
106 const poly integer_lcoeff () const;
107 const poly lcoeff_univar (int) const;
108 const poly lcoeff_multivar (int) const;
109 const poly coefficient (int, int) const;
110 const poly derivative (int) const;
111
112 // modulo calculus
113 void setmod(WORD, WORD=1);
114 void coefficients_modulo (UWORD *, WORD, bool);
115
116 // simple polynomials
117 static const poly simple_poly (PHEAD int, int=0, int=1, int=0, int=1);
118 static const poly simple_poly (PHEAD int, const poly&, int=1, int=0, int=1);
119
120 // conversion from/to form notation
121 static void get_variables (PHEAD std::vector<WORD *>, bool, bool);
122 static const poly argument_to_poly (PHEAD WORD *, bool, bool, poly *den=NULL);
123 static void poly_to_argument (const poly &, WORD *, bool);
124 static void poly_to_argument_with_den (const poly &, WORD, const UWORD *, WORD *, bool);
125 int size_of_form_notation () const;
126 int size_of_form_notation_with_den (WORD) const;
127 const poly & normalize ();
128
129 // operations for coefficient lists
130 static const poly from_coefficient_list (PHEAD const std::vector<WORD> &, int, WORD);
131 static const std::vector<WORD> to_coefficient_list (const poly &);
132 static const std::vector<WORD> coefficient_list_divmod (const std::vector<WORD> &, const std::vector<WORD> &, WORD, int);
133
134 // string output for debugging
135 const std::string to_string() const;
136
137 // monomials
138 static int monomial_compare (PHEAD const WORD *, const WORD *);
139 WORD last_monomial_index () const;
140 WORD* last_monomial () const;
141 int compare_degree_vector(const poly &) const;
142 std::vector<int> degree_vector() const;
143 int compare_degree_vector(const std::vector<int> &) const;
144
145 // (internal) mathematical operations
146 static void add (const poly &, const poly &, poly &);
147 static void sub (const poly &, const poly &, poly &);
148 static void mul (const poly &, const poly &, poly &);
149 static void div (const poly &, const poly &, poly &);
150 static void mod (const poly &, const poly &, poly &);
151 static void divmod (const poly &, const poly &, poly &, poly &, bool only_divides);
152 static bool divides (const poly &, const poly &);
153
154 static void mul_one_term (const poly &, const poly &, poly &);
155 static void mul_univar (const poly &, const poly &, poly &, int);
156 static void mul_heap (const poly &, const poly &, poly &);
157
158 static void divmod_one_term (const poly &, const poly &, poly &, poly &, bool);
159 static void divmod_univar (const poly &, const poly &, poly &, poly &, int, bool);
160 static void divmod_heap (const poly &, const poly &, poly &, poly &, bool, bool, bool&);
161
162 static void push_heap (PHEAD WORD **, int);
163 static void pop_heap (PHEAD WORD **, int);
164
165 PADPOINTER(1,0,2,0);
166 };
167
168 // comparison class for monomials (for std::sort)
169 class monomial_larger {
170
171 public:
172
173 #ifndef WITHPTHREADS
monomial_larger()174 monomial_larger() {}
175 #else
176 ALLPRIVATES *B;
177 monomial_larger(ALLPRIVATES *b): B(b) {}
178 #endif
179
operator()180 bool operator()(const WORD *a, const WORD *b) {
181 return poly::monomial_compare(BHEAD a, b) > 0;
182 }
183 };
184
185 // stream operator
186 std::ostream& operator<< (std::ostream &, const poly &);
187
188 // inline function definitions
189
190 /* Checks whether the terms array is large enough to add another
191 * term (of size AM.MaxTal) to the polynomials. In case not, it is
192 * expanded.
193 */
check_memory(int i)194 inline void poly::check_memory (int i) {
195 POLY_GETIDENTITY(*this);
196 if (i + 3 + AN.poly_num_vars + AM.MaxTal >= size_of_terms) expand_memory(i + AM.MaxTal);
197 // Used to be i+2 but there should also be space for a trailing zero
198 }
199
200 // indexing operators
201 inline WORD& poly::operator[] (int i) {
202 return terms[i];
203 }
204
205 inline const WORD& poly::operator[] (int i) const {
206 return terms[i];
207 }
208
209 /* Copies "num" WORD-sized terms from the pointer "source" to the
210 * current polynomial at index "dest"
211 */
termscopy(const WORD * source,int dest,int num)212 inline void poly::termscopy (const WORD *source, int dest, int num) {
213 memcpy (terms+dest, source, num*sizeof(WORD));
214 }
215