All Classes Functions Friends
mrpolynomial.h
1 #ifndef _MULTIPOLYNOMIAL_H_
2 #define _MULTIPOLYNOMIAL_H_
3 
4 /**
5  * Data Structure for multivariate rational number polynomial
6  * stored in a sparse case
7  **/
8 
9 
10 #include "../globals.h"
11 #include "../URPolynomial/urpolynomial.h"
12 #include "../Interval/interval.h"
13 #include "../../Polynomial/upolynomial.h"
14 
16  public:
17  int op; // 0: '+'; 1: '*'.
18  int type; // 0: coef & variate;
19  // 1: coef & result;
20  // 2: variate & result;
21  // 3: result & result;
22  int a; // variate or result or coef position
23  int b; // variate or result or coef position
24  Interval res; // Storing the result
25 };
26 
28  private:
29  int var; // Number of variables
30  std::string* names; // Variable names
31  // names[0] = "1" or "9" (user specified)
32  // names[1] < names[2] < ...
33  std::vector<RationalNumberTerm> terms; // Each term of the polynomial
34 
35  /* Return the first non-zero variate in the ascending order */
36  int firstNonZero(int* d);
37  /* Return the first non-equal variate in the descending order */
38  int checkDegs(int* d, int* b);
39  bool isOrderedRing(SparseMultivariateRationalPolynomial& b, std::vector<int>& xs);
40 
41  /* Real Root Isolation */
43  int refineSleeveUnivariateIntervals(Intervals*, Intervals*, Intervals*, DenseUnivariateRationalPolynomial*, DenseUnivariateRationalPolynomial*, lfixq);
44  void sleeveBoundURPolynomials(DenseUnivariateRationalPolynomial*, DenseUnivariateRationalPolynomial*, Intervals&, int, int);
45 
46  /**
47  * Compare two term degrees in the same Ring
48  *
49  * Return:
50  * 1: a > b
51  * 0: a = b
52  * -1: a < b
53  **/
54  int compareTermDegs(RationalNumberTerm& a, RationalNumberTerm& b);
55  int compareTermDegs(RationalNumberTerm& a, RationalNumberTerm& b, std::vector<int> xs);
56 
57  void basicOp(lfixq c, bool isOp);
58  void pomopo(RationalNumberTerm& t, SparseMultivariateRationalPolynomial& b, std::vector<int> xs);
59 
60  /* Is equal to another polynomial */
61  bool isEqual(SparseMultivariateRationalPolynomial& b);
62 
63  public:
64  static int characteristic;
65  static bool isPrimeField;
66  static bool isComplexField;
67  std::vector<SLPRepresentation> slp; // SLP's Representation
68 
69  /**
70  * Construct a multivariate polynomial
71  *
72  * @param
73  **/
75  names = new std::string[1];
76  names[0] = "1";
77  }
78  /**
79  * Construct a multivariate polynomial with number of terms and variates
80  *
81  * @param v: Number of variables
82  **/
84  var = v;
85  names = new std::string[var+1];
86  names[0] = "1";
87  for (int i = 1; i <= var; ++i) {
88  std::ostringstream convert;
89  convert << var - i + 1;
90  names[i] = "_";
91  names[i] += convert.str();
92  }
93  }
94  /**
95  * Construct with a variable name
96  * such that f(x) = x;
97  *
98  * @param x: The variable name
99  **/
101  var = 1;
102  names = new std::string[2];
103  names[0] = "9";
104  names[1] = x;
105  RationalNumberTerm t;
106  t.coef = 1;
107  t.v = 1;
108  t.degs = new int[1];
109  t.degs[0] = 1;
110  terms.push_back(t);
111  }
112  /**
113  * Copy Constructor
114  *
115  * @param b: A sparse multivariate polynomial
116  **/
117  SparseMultivariateRationalPolynomial(const SparseMultivariateRationalPolynomial& b) : var(b.var), terms(b.terms), slp(b.slp) {
118  names = new std::string[var+1];
119  std::copy(b.names, b.names+var+1, names);
120  }
121  /**
122  *
123  **/
125  var = 1;
126  names = new std::string[var+1];
127  names[0] = "9";
128  names[1] = p.variable();
129  for (int i = 0; i <= p.degree(); ++i) {
130  if (p.coefficient(i) != 0) {
131  RationalNumberTerm t;
132  t.coef = p.coefficient(i);
133  t.v = 1;
134  t.degs = new int[1];
135  t.degs[0] = i;
136  terms.push_back(t);
137  }
138  }
139  }
140  /**
141  * Construct from a SUP<SMQP> polynomial
142  *
143  * @param s: The SUP<SMQP> polynomial
144  **/
146  /**
147  * Destroy the polynomial
148  *
149  * @param
150  **/
152  slp.clear();
153  terms.clear();
154  delete [] names;
155  }
156 
157  /**
158  * Get number of variables
159  *
160  * @param
161  **/
162  inline int numberOfVariables() {
163  return var;
164  }
165 
166  /**
167  * Get number of non-zero terms
168  *
169  * @param
170  **/
171  inline int numberOfTerms() {
172  return terms.size();
173  }
174  /**
175  * Get the degree of a variable
176  *
177  * @param x: The variable name
178  **/
179  inline int degree(std::string x) {
180  int k = 0;
181  for (int i = 1; i <= var; ++i) {
182  if (names[i] == x) {
183  k = i;
184  break;
185  }
186  }
187  if (!k) { return 0; }
188  int d = 0;
189  for (int i = terms.size()-1; i > -1; --i) {
190  if (terms[i].degs[k] > d)
191  d = terms[i].degs[k];
192  }
193  return d;
194  }
195  /**
196  * Get the leading variable
197  *
198  * @param
199  **/
200  inline std::string leadingVariable() {
201  return names[var];
202  }
203  /**
204  * Get the leading term's leading variable's degree
205  *
206  * @param
207  **/
208  inline int leadingVariableDegree() {
209  int n = terms.size();
210  if (n)
211  return terms[n-1].degs[var-1];
212  else
213  return 0;
214  }
215  /**
216  * Get the leading coefficient
217  *
218  * @param
219  **/
220  inline mpq_class leadingCoefficient() {
221  int n = terms.size();
222  if (n)
223  return terms[n].coef;
224  else
225  return 0;
226  }
227  /**
228  * Get the leading coefficient over a variable
229  *
230  * @param x: The name of the variable
231  * @param e: The leading exponent of the variable
232  **/
234  /**
235  * Is zero polynomial
236  *
237  * @param
238  **/
239  inline bool isZero() {
240  return !terms.size();
241  }
242  /**
243  * Zero polynomial
244  *
245  * @param
246  **/
247  inline void zero() {
248  terms.clear();
249  }
250  /**
251  * Is polynomial a constant 1
252  *
253  * @param
254  **/
255  inline bool isOne() {
256  if (terms.size() == 1) {
257  for (int i = 0; i < var; ++i) {
258  if (terms[0].degs[i] != 0)
259  return 0;
260  }
261  if (terms[0].coef == 1)
262  return 1;
263  }
264  return 0;
265  }
266  /**
267  * Set polynomial to 1
268  *
269  * @param
270  **/
271  inline void one() {
272  terms.clear();
273  RationalNumberTerm t;
274  t.coef = 1;
275  t.v = var;
276  t.degs = new int[var];
277  for (int i = 0; i < var; ++i)
278  t.degs[i] = 0;
279  terms.push_back(t);
280  }
281  /**
282  * Is polynomial a constant -1
283  *
284  * @param
285  **/
286  inline bool isNegativeOne() {
287  if (terms.size() == 1) {
288  for (int i = 0; i < var; ++i) {
289  if (terms[0].degs[i] != 0)
290  return 0;
291  }
292  if (terms[0].coef == -1)
293  return 1;
294  }
295  return 0;
296  }
297  /**
298  * Set polynomial to -1
299  *
300  * @param
301  **/
302  inline void negativeOne() {
303  terms.clear();
304  RationalNumberTerm t;
305  t.coef = -1;
306  t.v = var;
307  t.degs = new int[var];
308  for (int i = 0; i < var; ++i)
309  t.degs[i] = 0;
310  terms.push_back(t);
311  }
312  /**
313  * Is a constant
314  *
315  * @param
316  **/
317  inline int isConstant() {
318  if (!terms.size())
319  return 1;
320  else if (terms.size() == 1) {
321  for (int i = 0; i < var; ++i) {
322  if (terms[0].degs[i] != 0)
323  return 0;
324  }
325  if (terms[0].coef >= 0) { return 1; }
326  else { return -1; }
327  }
328  return 0;
329  }
330  /**
331  * Overload operator =
332  *
333  * @param b: A sparse multivariate polynomial
334  **/
336  if (this != &b) {
337  terms.clear();
338  delete [] names;
339 
340  var = b.numberOfVariables();
341  names = new std::string[var+1];
342  std::copy(b.names, b.names+var+1, names);
343  terms = b.terms;
344  slp = b.slp;
345  }
346  return *this;
347  }
348  /**
349  * Overload operator ==
350  *
351  * @param b: A multivariate rational polynomial
352  **/
354  return isEqual(b);
355  }
356  /**
357  * Overload operator !=
358  *
359  * @param b: A multivariate rational polynomial
360  **/
362  return !isEqual(b);
363  }
364  /**
365  * Overload operator ^
366  * replace xor operation by exponentiation
367  *
368  * @param e: The exponentiation, e > 0
369  **/
371  /**
372  * Overload operator ^=
373  * replace xor operation by exponentiation
374  *
375  * @param e: The exponentiation, e > 0
376  **/
378  *this = *this ^ e;
379  return *this;
380  }
381 
382  /**
383  * Overload operator +
384  *
385  * @param b: A multivariate rational polynomial
386  **/
388  /**
389  * Overload operator +=
390  *
391  * @param b: A multivariate rational polynomial
392  **/
394  *this = *this + b;
395  return *this;
396  }
397  /**
398  * Overload operator +
399  *
400  * @param c: A constant
401  **/
404  return (r += c);
405  }
408  return (r += c);
409  }
410  /**
411  * Overload operator +=
412  *
413  * @param c: A constant
414  **/
416  basicOp(c, 0);
417  return *this;
418  }
420  mpq_class e = lfixq(c.get_mpq_t());
421  return (*this += e);
422  }
424  return (p + c);
425  }
426  /**
427  * Overload operator -
428  *
429  * @param b: A multivariate rational polynomial
430  **/
432  /**
433  * Overload operator -=
434  *
435  * @param b: A multivariate rational polynomial
436  **/
438  *this = *this - b;
439  return *this;
440  }
441  /**
442  * Overload operator -
443  *
444  * @param c: A constant
445  **/
448  return (r -= c);
449  }
452  return (r -= c);
453  }
454  /**
455  * Overload operator -=
456  *
457  * @param c: A constant
458  **/
460  basicOp(c, 1);
461  return *this;
462  }
464  mpq_class e = lfixq(c.get_mpq_t());
465  return (*this += e);
466  }
467 
469  return -p + c;
470  }
471  /**
472  * Overload operator - (negate)
473  *
474  * @param
475  **/
477  /**
478  * Overload operator *
479  *
480  * @param b: A multivariate rational polynomial
481  **/
483  /**
484  * Overload operator *=
485  *
486  * @param b: A multivariate rational polynomial
487  **/
489  *this = *this * b;
490  return *this;
491  }
492  /**
493  * Overload operator *
494  *
495  * @param e: A constant
496  **/
499  return (r *= e);
500  }
503  return (r *= e);
504  }
507  return (r *= e);
508  }
509  /**
510  * Overload operator *=
511  *
512  * @param e: A constant
513  **/
515  if (e != 0 && e != 1) {
516  for (int i = 0; i < terms.size(); ++i)
517  terms[i].coef *= e;
518  }
519  else if (e == 0) { terms.clear(); }
520  return *this;
521  }
523  mpq_class c = lfixq(e.get_mpq_t());
524  return (*this += c);
525  }
526 
528  return (p * c);
529  }
531  return (p * c);
532  }
533  /**
534  * Overload operator *=
535  *
536  * @param e: A machine-word constant
537  **/
539  if (e != 0 && e != 1) {
540  for (int i = 0; i < terms.size(); ++i)
541  terms[i].coef *= e;
542  }
543  else if (e == 0) { terms.clear(); }
544  return *this;
545  }
546  /**
547  * Overload operator /
548  *
549  * @param b: A multivariate polynomial
550  **/
553  return (r /= b);
554  }
555  /**
556  * Overload operator /=
557  *
558  * @param b: A multivariate polynomial
559  **/
561  /**
562  * Overload operator /
563  *
564  * @param e: A constant
565  **/
568  return (r /= e);
569  }
572  return (r /= e);
573  }
574  /**
575  * Overload operator /=
576  *
577  * @param e: A constant
578  **/
580  if (e == 0) {
581  std::cout << "BPAS: error, dividend is zero." << std::endl;
582  exit(1);
583  }
584  else if (e != 1) {
585  for (int i = 0; i < terms.size(); ++i)
586  terms[i].coef /= e;
587  }
588  return *this;
589  }
591  mpq_class c = lfixq(e.get_mpq_t());
592  return (*this += c);
593  }
594 
596  /**
597  * Set variable names
598  *
599  * @param xs: Variable names
600  **/
601  void setVariableNames (std::vector<std::string> xs);
602  /**
603  * Get variable names
604  *
605  * @param
606  **/
607  inline std::vector<std::string> variables() {
608  std::vector<std::string> xs;
609  for (int i = var; i > 0; --i)
610  xs.push_back(names[i]);
611  return xs;
612  }
613  /**
614  * Get the coefficient of one term
615  *
616  * @param d: The exponet of each variable
617  * @param v: Number of variables
618  **/
619  mpq_class coefficient(int v, int* d);
620  /**
621  * Set the coefficients of one term
622  *
623  * @param v: Number of variables
624  * @param d: Its exponent of each variable
625  * @param val: Coefficient
626  **/
627  void setCoefficient(int v, int* d, mpq_class val);
628 
629  /**
630  * Negate all the coefficients
631  *
632  * @param
633  **/
634  void negate();
635 
636  /**
637  * Has the zero root of the most weighted variate
638  *
639  * @param
640  **/
642  /**
643  * Convert to SUP<SMQP>
644  *
645  * @param x: The variable name
646  **/
648 
649  /**
650  * SLP's representation of the polynomial
651  *
652  * @param
653  **/
654  void straightLineProgram();
655  /**
656  * Print SLP representation
657  *
658  * @param
659  **/
660  void printSLP();
661  /**
662  * Given one real root for x_1, .., x_{n-1},
663  * isolate positive roots for x_n
664  *
665  * @param mpIs: Roots of x_n (Output)
666  * @param apIs: A root of previous polynomials
667  * @param s: deal with 0: positive roots; 1: negative roots
668  * @param check: 1: check the leading or tail coefficient; 0: do not
669  * @param width: Interval's right - left < width
670  * @param ts: Taylor Shift option
671  *
672  * Return
673  * 1: Need to refine preious polynomials
674  * 0: Found positive real roots
675  **/
676  int positiveRealRootIsolation(Intervals* pIs, Intervals& apIs, mpq_class width, int ts=-1, bool s=0, bool check=1);
677  /**
678  * Overload stream operator <<
679  *
680  * @param out: Stream object
681  * @param b: A multivariate rational polynoial
682  **/
683  friend std::ostream& operator<< (std::ostream& out, SparseMultivariateRationalPolynomial b);
684 };
685 
686 #endif
687 /* This file is part of the BPAS library http://www.bpaslib.org
688 
689  BPAS is free software: you can redistribute it and/or modify
690  it under the terms of the GNU General Public License as published by
691  the Free Software Foundation, either version 3 of the License, or
692  (at your option) any later version.
693 
694  BPAS is distributed in the hope that it will be useful,
695  but WITHOUT ANY WARRANTY; without even the implied warranty of
696  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
697  GNU General Public License for more details.
698 
699  You should have received a copy of the GNU General Public License
700  along with BPAS. If not, see <http://www.gnu.org/licenses/>.
701 
702  Copyright: Changbo Chen <changbo.chen@hotmail.com>
703  Farnam Mansouri <mansouri.farnam@gmail.com>
704  Marc Moreno Maza <moreno@csd.uwo.ca>
705  Ning Xie <nxie6@csd.uwo.ca>
706  Yuzhen Xie <yuzhenxie@yahoo.ca>
707 
708 */
709 
710 
int numberOfVariables()
Definition: mrpolynomial.h:162
int positiveRealRootIsolation(Intervals *pIs, Intervals &apIs, mpq_class width, int ts=-1, bool s=0, bool check=1)
Definition: upolynomial.h:10
std::string variable()
Definition: urpolynomial.h:206
Definition: mrpolynomial.h:27
SparseMultivariateRationalPolynomial operator*(SparseMultivariateRationalPolynomial b)
mpq_class leadingCoefficient()
Definition: mrpolynomial.h:220
bool operator!=(SparseMultivariateRationalPolynomial &b)
Definition: mrpolynomial.h:361
SparseMultivariateRationalPolynomial & operator-=(SparseMultivariateRationalPolynomial b)
Definition: mrpolynomial.h:437
mpq_class coefficient(int v, int *d)
void negativeOne()
Definition: mrpolynomial.h:302
Definition: interval.h:14
SparseMultivariateRationalPolynomial(int v)
Definition: mrpolynomial.h:83
SparseMultivariateRationalPolynomial(std::string x)
Definition: mrpolynomial.h:100
int leadingVariableDegree()
Definition: mrpolynomial.h:208
SparseUnivariatePolynomial< SparseMultivariateRationalPolynomial > convertToSUP(std::string x)
Definition: urpolynomial.h:16
bool isZero()
Definition: mrpolynomial.h:239
Definition: polynomial.h:74
friend std::ostream & operator<<(std::ostream &out, SparseMultivariateRationalPolynomial b)
int isConstant()
Definition: mrpolynomial.h:317
SparseMultivariateRationalPolynomial & operator+=(SparseMultivariateRationalPolynomial b)
Definition: mrpolynomial.h:393
SparseMultivariateRationalPolynomial()
Definition: mrpolynomial.h:74
void zero()
Definition: mrpolynomial.h:247
SparseMultivariateRationalPolynomial operator/(SparseMultivariateRationalPolynomial b)
Definition: mrpolynomial.h:551
~SparseMultivariateRationalPolynomial()
Definition: mrpolynomial.h:151
void setCoefficient(int v, int *d, mpq_class val)
SparseMultivariateRationalPolynomial & operator*=(SparseMultivariateRationalPolynomial b)
Definition: mrpolynomial.h:488
bool operator==(SparseMultivariateRationalPolynomial &b)
Definition: mrpolynomial.h:353
std::vector< std::string > variables()
Definition: mrpolynomial.h:607
SparseMultivariateRationalPolynomial & operator^=(int e)
Definition: mrpolynomial.h:377
int degree(std::string x)
Definition: mrpolynomial.h:179
void one()
Definition: mrpolynomial.h:271
SparseMultivariateRationalPolynomial & operator/=(SparseMultivariateRationalPolynomial b)
bool isNegativeOne()
Definition: mrpolynomial.h:286
Definition: ring.h:166
SparseMultivariateRationalPolynomial operator+(SparseMultivariateRationalPolynomial b)
Definition: mrpolynomial.h:15
SparseMultivariateRationalPolynomial leadingCoefficientInVariable(std::string x, int *e=NULL)
mpq_class coefficient(int k)
Definition: urpolynomial.h:174
SparseMultivariateRationalPolynomial(const SparseMultivariateRationalPolynomial &b)
Definition: mrpolynomial.h:117
bool isOne()
Definition: mrpolynomial.h:255
std::string leadingVariable()
Definition: mrpolynomial.h:200
Definition: interval.h:35
SparseMultivariateRationalPolynomial & operator=(SparseMultivariateRationalPolynomial b)
Definition: mrpolynomial.h:335
SparseMultivariateRationalPolynomial operator^(int e)
int numberOfTerms()
Definition: mrpolynomial.h:171
void setVariableNames(std::vector< std::string > xs)
SparseMultivariateRationalPolynomial operator-()
int degree()
Definition: urpolynomial.h:145