3 * Some test with polynomial GCD calculations. See also the checks for
4 * rational function normalization in normalization.cpp. */
7 * GiNaC Copyright (C) 1999 Johannes Gutenberg University Mainz, Germany
9 * This program is free software; you can redistribute it and/or modify
10 * it under the terms of the GNU General Public License as published by
11 * the Free Software Foundation; either version 2 of the License, or
12 * (at your option) any later version.
14 * This program is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 * GNU General Public License for more details.
19 * You should have received a copy of the GNU General Public License
20 * along with this program; if not, write to the Free Software
21 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
24 #include <ginac/ginac.h>
26 #ifndef NO_GINAC_NAMESPACE
27 using namespace GiNaC;
28 #endif // ndef NO_GINAC_NAMESPACE
30 const int MAX_VARIABLES = 5;
32 static symbol x("x"), z("z");
33 static symbol y[MAX_VARIABLES];
36 static unsigned poly_gcd1(void)
38 for (int v=1; v<=MAX_VARIABLES; v++) {
41 for (int i=0; i<v; i++) {
46 ex f = (e1 + 1) * (e1 + 2);
47 ex g = e2 * (-pow(x, 2) * y[0] * 3 + pow(y[0], 2) - 1);
50 clog << "case 1, gcd(" << f << "," << g << ") = " << r << " (should be 1)" << endl;
57 // Linearly dense quartic inputs with quadratic GCDs
58 static unsigned poly_gcd2(void)
60 for (int v=1; v<=MAX_VARIABLES; v++) {
63 for (int i=0; i<v; i++) {
68 ex d = pow(e1 + 1, 2);
69 ex f = d * pow(e2 - 2, 2);
70 ex g = d * pow(e1 + 2, 2);
75 if ((r - d).expand().compare(exZERO()) != 0) {
76 clog << "case 2, gcd(" << f << "," << g << ") = " << r << " (should be " << d << ")" << endl;
83 // Sparse GCD and inputs where degrees are proportional to the number of variables
84 static unsigned poly_gcd3(void)
86 for (int v=1; v<=MAX_VARIABLES; v++) {
87 ex e1 = pow(x, v + 1);
88 for (int i=0; i<v; i++)
89 e1 += pow(y[i], v + 1);
95 if ((r - d).expand().compare(exZERO()) != 0) {
96 clog << "case 3, gcd(" << f << "," << g << ") = " << r << " (should be " << d << ")" << endl;
103 // Variation of case 3; major performance degradation with PRS
104 static unsigned poly_gcd3p(void)
106 for (int v=1; v<=MAX_VARIABLES; v++) {
107 ex e1 = pow(x, v + 1);
109 for (int i=0; i<v; i++) {
110 e1 += pow(y[i], v + 1);
118 if ((r - d).expand().compare(exZERO()) != 0) {
119 clog << "case 3p, gcd(" << f << "," << g << ") = " << r << " (should be " << d << ")" << endl;
126 // Quadratic non-monic GCD; f and g have other quadratic factors
127 static unsigned poly_gcd4(void)
129 for (int v=1; v<=MAX_VARIABLES; v++) {
130 ex e1 = pow(x, 2) * pow(y[0], 2);
131 ex e2 = pow(x, 2) - pow(y[0], 2);
133 for (int i=1; i<v; i++) {
141 ex g = d * pow(e3 + 2, 2);
143 if ((r - d).expand().compare(exZERO()) != 0) {
144 clog << "case 4, gcd(" << f << "," << g << ") = " << r << " (should be " << d << ")" << endl;
151 // Completely dense non-monic quadratic inputs with dense non-monic linear GCDs
152 static unsigned poly_gcd5(void)
154 for (int v=1; v<=MAX_VARIABLES; v++) {
158 for (int i=0; i<v; i++) {
168 if ((r - d).expand().compare(exZERO()) != 0) {
169 clog << "case 5, gcd(" << f << "," << g << ") = " << r << " (should be " << d << ")" << endl;
176 // Sparse non-monic quadratic inputs with linear GCDs
177 static unsigned poly_gcd5p(void)
179 for (int v=1; v<=MAX_VARIABLES; v++) {
181 for (int i=0; i<v; i++)
188 if ((r - d).expand().compare(exZERO()) != 0) {
189 clog << "case 5p, gcd(" << f << "," << g << ") = " << r << " (should be " << d << ")" << endl;
196 // Trivariate inputs with increasing degrees
197 static unsigned poly_gcd6(void)
201 for (int j=1; j<=MAX_VARIABLES; j++) {
202 ex d = pow(x, j) * y * (z - 1);
203 ex f = d * (pow(x, j) + pow(y, j + 1) * pow(z, j) + 1);
204 ex g = d * (pow(x, j + 1) + pow(y, j) * pow(z, j + 1) - 7);
206 if ((r - d).expand().compare(exZERO()) != 0) {
207 clog << "case 6, gcd(" << f << "," << g << ") = " << r << " (should be " << d << ")" << endl;
214 // Trivariate polynomials whose GCD has common factors with its cofactors
215 static unsigned poly_gcd7(void)
218 ex p = x - y * z + 1;
219 ex q = x - y + z * 3;
221 for (int j=1; j<=3; j++) {
222 for (int k=j+1; k<=4; k++) {
223 ex d = pow(p, j) * pow(q, j);
224 ex f = pow(p, j) * pow(q, k);
225 ex g = pow(p, k) * pow(q, j);
227 if ((r - d).expand().compare(exZERO()) != 0 && (r + d).expand().compare(exZERO()) != 0) {
228 clog << "case 7, gcd(" << f << "," << g << ") = " << r << " (should be " << d << ")" << endl;
236 unsigned poly_gcd(void)
240 cout << "checking polynomial GCD computation..." << flush;
241 clog << "---------polynomial GCD computation:" << endl;
243 result += poly_gcd1();
244 result += poly_gcd2();
245 result += poly_gcd3();
246 // result += poly_gcd3p(); // takes extremely long (PRS "worst" case)
247 result += poly_gcd4();
248 result += poly_gcd5();
249 result += poly_gcd5p();
250 result += poly_gcd6();
251 result += poly_gcd7();
255 clog << "(no output)" << endl;