source: git/libpolys/coeffs/flintcf_Qrat.cc @ a507140

spielwiese
Last change on this file since a507140 was a507140, checked in by Hans Schoenemann <hannes@…>, 4 years ago
opt: nCoeffWrite/nCoeffString via nCoeffName
  • Property mode set to 100644
File size: 44.8 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5* ABSTRACT: flint: rational functions over Q (using fmpq_mpoly)
6*/
7#include <ctype.h> /* isdigit*/
8
9#include "misc/auxiliary.h"
10
11#include "coeffs/coeffs.h"
12
13#ifdef HAVE_FLINT
14#include "flint/flint.h"
15#if __FLINT_RELEASE >= 20503
16#include "factory/factory.h"
17
18#include "coeffs/numbers.h"
19#include "coeffs/longrat.h"
20#include "coeffs/flintcf_Qrat.h"
21#include "polys/flint_mpoly.h"
22#ifdef QA_DEBUG
23#define TRANSEXT_PRIVATES
24#include "polys/ext_fields/transext.h"
25#include "polys/monomials/p_polys.h"
26#endif
27
28typedef fmpq_rat_struct *fmpq_rat_ptr;
29typedef fmpq_mpoly_struct *fmpq_mpoly_ptr;
30typedef fmpq_mpoly_ctx_struct *fmpq_ctx_ptr;
31typedef fmpz *fmpz_ptr;
32typedef fmpq_rat_data_struct *data_ptr;
33
34/******************************************************************************
35* Helper functions
36******************************************************************************/
37
38/*2
39* extracts a long integer from s, returns the rest
40*/
41static char * nlEatLong(char *s, fmpz_ptr i)
42{
43  const char * start = s;
44
45  while (*s >= '0' && *s <= '9') s++;
46  if (*s == '\0')
47  {
48    fmpz_set_str(i, start, 10);
49  }
50  else
51  {
52    char c = *s;
53    *s = '\0';
54    fmpz_set_str(i, start, 10);
55    *s = c;
56  }
57  return s;
58}
59
60static void fmpq_rat_init(fmpq_rat_ptr a, const coeffs r)
61{
62  fmpq_mpoly_init(a->num, ((data_ptr)r->data)->ctx);
63  fmpq_mpoly_init(a->den, ((data_ptr)r->data)->ctx);
64}
65
66static void fmpq_rat_clear(fmpq_rat_ptr a, const coeffs r)
67{
68  fmpq_mpoly_clear(a->num, ((data_ptr)r->data)->ctx);
69  fmpq_mpoly_clear(a->den, ((data_ptr)r->data)->ctx);
70}
71
72static void fmpq_rat_canonicalise(fmpq_rat_ptr a, const coeffs r)
73{
74  fmpz_t n, d;
75  fmpz_init(n);
76  fmpz_init(d);
77  fmpz_gcd(n, fmpq_numref(a->num->content), fmpq_numref(a->den->content));
78  fmpz_lcm(d, fmpq_denref(a->num->content), fmpq_denref(a->den->content));
79  if (!fmpz_is_one(d))
80  {
81     fmpq_mul_fmpz(a->num->content, a->num->content, d);
82     fmpq_mul_fmpz(a->den->content, a->den->content, d);
83  }
84  if (!fmpz_is_one(n))
85  {
86     fmpq_div_fmpz(a->num->content, a->num->content, n);
87     fmpq_div_fmpz(a->den->content, a->den->content, n);
88  }
89  fmpz_clear(n);
90  fmpz_clear(d);
91}
92
93/******************************************************************************
94* Main interface
95******************************************************************************/
96
97static BOOLEAN CoeffIsEqual(const coeffs c, n_coeffType n, void * parameter)
98{
99  if (c->type == n)
100  {
101    const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
102    const QaInfo *par=(QaInfo*)parameter;
103    if (par->N != c->iNumberOfParameters) return FALSE;
104    // compare parameter names
105    for(int i=0;i<par->N;i++)
106    {
107      if (strcmp(par->names[i],c->pParameterNames[i])!=0) return FALSE;
108    }
109    return TRUE;
110  }
111  return FALSE;
112}
113
114static number Mult(number a, number b, const coeffs c)
115{
116  n_Test(a,c);
117  n_Test(b,c);
118  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
119  fmpq_rat_init(res, c);
120  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
121  const fmpq_rat_ptr y = (fmpq_rat_ptr) b;
122  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
123  if (fmpq_mpoly_equal(x->den, y->den, ctx)) /* denominators equal */
124  {
125    fmpq_mpoly_mul(res->num, x->num, y->num, ctx);
126    fmpq_mpoly_mul(res->den, x->den, y->den, ctx);
127  }
128  else if (fmpq_mpoly_is_one(x->den, ctx)) /* first denominator 1 */
129  {
130    fmpq_mpoly_t gd;
131    fmpq_mpoly_init(gd, ctx);
132    fmpq_mpoly_gcd(gd, x->num, y->den, ctx);
133    if (fmpq_mpoly_is_one(gd, ctx))
134    {
135      fmpq_mpoly_mul(res->num, x->num, y->num, ctx);
136      fmpq_mpoly_set(res->den, y->den, ctx);
137    }
138    else
139    {
140      fmpq_mpoly_div(res->num, x->num, gd, ctx);
141      fmpq_mpoly_mul(res->num, res->num, y->num, ctx);
142      fmpq_mpoly_div(res->den, y->den, gd, ctx);
143    }
144    fmpq_mpoly_clear(gd, ctx);
145  }
146  else if (fmpq_mpoly_is_one(y->den, ctx)) /* second denominator 1 */
147  {
148    fmpq_mpoly_t gd;
149    fmpq_mpoly_init(gd, ctx);
150    fmpq_mpoly_gcd(gd, y->num, x->den, ctx);
151    if (fmpq_mpoly_is_one(gd, ctx))
152    {
153      fmpq_mpoly_mul(res->num, x->num, y->num, ctx);
154      fmpq_mpoly_set(res->den, x->den, ctx);
155    }
156    else
157    {
158      fmpq_mpoly_div(res->num, y->num, gd, ctx);
159      fmpq_mpoly_mul(res->num, res->num, x->num, ctx);
160      fmpq_mpoly_div(res->den, x->den, gd, ctx);
161    }
162    fmpq_mpoly_clear(gd, ctx);
163  }
164  else /* general case */
165  {
166    fmpq_mpoly_t g1, g2;
167    fmpq_mpoly_ptr n1, n2, d1, d2;
168    fmpq_mpoly_init(g1, ctx);
169    fmpq_mpoly_init(g2, ctx);
170    fmpq_mpoly_gcd(g1, x->num, y->den, ctx);
171    fmpq_mpoly_gcd(g2, y->num, x->den, ctx);
172    n1 = x->num; d2 = y->den;
173    d1 = x->den; n2 = y->num;
174    if (!fmpq_mpoly_is_one(g1, ctx))
175    {
176      fmpq_mpoly_div(res->num, x->num, g1, ctx);
177      fmpq_mpoly_div(g1, y->den, g1, ctx);
178      n1 = res->num; d2 = g1;
179    }
180    if (!fmpq_mpoly_is_one(g2, ctx))
181    {
182      fmpq_mpoly_div(res->den, y->num, g2, ctx);
183      fmpq_mpoly_div(g2, x->den, g2, ctx);
184      n2 = res->den; d1 = g2;
185    }
186    fmpq_mpoly_mul(res->num, n1, n2, ctx);
187    fmpq_mpoly_mul(res->den, d1, d2, ctx);
188    fmpq_mpoly_clear(g1, ctx);
189    fmpq_mpoly_clear(g2, ctx);
190  }
191  fmpq_rat_canonicalise(res, c);
192  #ifdef QA_DEBUG
193  res->p=n_Mult(x->p,y->p, ((data_ptr)c->data)->C);
194  #endif
195  n_Test((number)res, c);
196  return (number) res;
197}
198
199static number Sub(number a, number b, const coeffs c)
200{
201  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
202  fmpq_rat_init(res, c);
203  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
204  const fmpq_rat_ptr y = (fmpq_rat_ptr) b;
205  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
206  if (fmpq_mpoly_equal(x->den, y->den, ctx)) /* denominators equal */
207  {
208    fmpq_mpoly_sub(res->num, x->num, y->num, ctx);
209    if (fmpq_mpoly_is_zero(res->num, ctx))
210    {
211      fmpq_mpoly_one(res->den, ctx);
212      n_Test((number)res,c);
213      return (number)res;
214    }
215    else
216    if (fmpq_mpoly_is_one(x->den, ctx))
217    {
218      fmpq_mpoly_set(res->den, x->den, ctx);
219      n_Test((number)res,c);
220      return (number)res;
221    }
222    else
223    {
224      fmpq_mpoly_t gd;
225      fmpq_mpoly_init(gd, ctx);
226      fmpq_mpoly_gcd(gd, res->num, x->den, ctx);
227      if (fmpq_mpoly_is_one(gd, ctx))
228      {
229        fmpq_mpoly_set(res->den, x->den, ctx);
230      }
231      else
232      {
233        fmpq_mpoly_div(res->den, x->den, gd, ctx);
234        fmpq_mpoly_div(res->num, res->num, gd, ctx);
235      }
236      fmpq_mpoly_clear(gd, ctx);
237    }
238  }
239  else if (fmpq_mpoly_is_one(x->den, ctx)) /* first denominator 1 */
240  {
241    fmpq_mpoly_mul(res->num, x->num, y->den, ctx);
242    fmpq_mpoly_sub(res->num, res->num, y->num, ctx);
243    if (fmpq_mpoly_is_zero(res->num, ctx))
244    {
245      fmpq_mpoly_one(res->den, ctx);
246      n_Test((number)res,c);
247      return (number)res;
248    }
249    else
250    {
251      fmpq_mpoly_set(res->den, y->den, ctx);
252    }
253  }
254  else if (fmpq_mpoly_is_one(y->den, ctx)) /* second denominator 1 */
255  {
256    fmpq_mpoly_mul(res->num, y->num, x->den, ctx);
257    fmpq_mpoly_sub(res->num, x->num, res->num, ctx);
258    if (fmpq_mpoly_is_zero(res->num,ctx))
259    {
260      fmpq_mpoly_one(res->den, ctx);
261      n_Test((number)res,c);
262      return (number)res;
263    }
264    else
265    {
266      fmpq_mpoly_set(res->den, x->den, ctx);
267    }
268  }
269  else /* general case */
270  {
271    fmpq_mpoly_t gd;
272    fmpq_mpoly_init(gd, ctx);
273    fmpq_mpoly_gcd(gd, x->den, y->den, ctx);
274    if (fmpq_mpoly_is_one(gd, ctx))
275    {
276      fmpq_mpoly_mul(res->num, x->num, y->den, ctx);
277      fmpq_mpoly_mul(gd, y->num, x->den, ctx);
278      fmpq_mpoly_sub(res->num, res->num, gd, ctx);
279      if (fmpq_mpoly_is_zero(res->num,ctx))
280      {
281        fmpq_mpoly_one(res->den, ctx);
282        n_Test((number)res,c);
283        return (number)res;
284      }
285      else
286      {
287        fmpq_mpoly_mul(res->den, x->den, y->den, ctx);
288      }
289    }
290    else
291    {
292      fmpq_mpoly_t q2;
293      fmpq_mpoly_init(q2, ctx);
294      fmpq_mpoly_div(res->den, x->den, gd, ctx);
295      fmpq_mpoly_div(q2, y->den, gd, ctx);
296      fmpq_mpoly_mul(res->num, q2, x->num, ctx);
297      fmpq_mpoly_mul(res->den, res->den, y->num, ctx);
298      fmpq_mpoly_sub(res->num, res->num, res->den, ctx);
299      fmpq_mpoly_gcd(res->den, res->num, gd, ctx);
300      if (fmpq_mpoly_is_one(res->den, ctx))
301      {
302        fmpq_mpoly_mul(res->den, q2, x->den, ctx);
303      }
304      else
305      {
306        fmpq_mpoly_div(res->num, res->num, res->den, ctx);
307        fmpq_mpoly_div(gd, x->den, res->den, ctx);
308        fmpq_mpoly_mul(res->den, gd, q2, ctx);
309      }
310      fmpq_mpoly_clear(q2, ctx);
311    }
312    fmpq_mpoly_clear(gd, ctx);
313  }
314  #ifdef QA_DEBUG
315  res->p=n_Sub(x->p,y->p, ((data_ptr)c->data)->C);
316  #endif
317  n_Test((number)res, c);
318  return (number) res;
319}
320
321static number Add(number a, number b, const coeffs c)
322{
323  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
324  fmpq_rat_init(res, c);
325  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
326  const fmpq_rat_ptr y = (fmpq_rat_ptr) b;
327  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
328  if (fmpq_mpoly_equal(x->den, y->den, ctx)) /* denominators equal */
329  {
330    fmpq_mpoly_add(res->num, x->num, y->num, ctx);
331    if (fmpq_mpoly_is_zero(res->num, ctx))
332    {
333      fmpq_mpoly_one(res->den, ctx);
334      n_Test((number)res,c);
335      return (number)res;
336    }
337    else
338    if (fmpq_mpoly_is_one(x->den, ctx))
339    {
340      fmpq_mpoly_set(res->den, x->den, ctx);
341      n_Test((number)res,c);
342      return (number)res;
343    }
344    else
345    {
346      fmpq_mpoly_t gd;
347      fmpq_mpoly_init(gd, ctx);
348      fmpq_mpoly_gcd(gd, res->num, x->den, ctx);
349      if (fmpq_mpoly_is_one(gd, ctx))
350      {
351        fmpq_mpoly_set(res->den, x->den, ctx);
352      }
353      else
354      {
355        fmpq_mpoly_div(res->den, x->den, gd, ctx);
356        fmpq_mpoly_div(res->num, res->num, gd, ctx);
357      }
358      fmpq_mpoly_clear(gd, ctx);
359    }
360  }
361  else if (fmpq_mpoly_is_one(x->den, ctx)) /* first denominator 1 */
362  {
363    fmpq_mpoly_mul(res->num, x->num, y->den, ctx);
364    fmpq_mpoly_add(res->num, res->num, y->num, ctx);
365    if (fmpq_mpoly_is_zero(res->num, ctx))
366    {
367      fmpq_mpoly_one(res->den, ctx);
368      n_Test((number)res,c);
369      return (number)res;
370    }
371    else
372    {
373      fmpq_mpoly_set(res->den, y->den, ctx);
374    }
375  }
376  else if (fmpq_mpoly_is_one(y->den, ctx)) /* second denominator 1 */
377  {
378    fmpq_mpoly_mul(res->num, y->num, x->den, ctx);
379    fmpq_mpoly_add(res->num, x->num, res->num, ctx);
380    if (fmpq_mpoly_is_zero(res->num, ctx))
381    {
382      fmpq_mpoly_one(res->den, ctx);
383      n_Test((number)res,c);
384      return (number)res;
385    }
386    else
387    {
388      fmpq_mpoly_set(res->den, x->den, ctx);
389    }
390  }
391  else /* general case */
392  {
393    fmpq_mpoly_t gd;
394    fmpq_mpoly_init(gd, ctx);
395    fmpq_mpoly_gcd(gd, x->den, y->den, ctx);
396    if (fmpq_mpoly_is_one(gd, ctx))
397    {
398      fmpq_mpoly_mul(res->num, x->num, y->den, ctx);
399      fmpq_mpoly_mul(gd, y->num, x->den, ctx);
400      fmpq_mpoly_add(res->num, res->num, gd, ctx);
401      if (fmpq_mpoly_is_zero(res->num,ctx))
402      {
403        fmpq_mpoly_one(res->den, ctx);
404        n_Test((number)res,c);
405        return (number)res;
406      }
407      else
408      {
409        fmpq_mpoly_mul(res->den, x->den, y->den, ctx);
410      }
411    }
412    else
413    {
414      fmpq_mpoly_t q2;
415      fmpq_mpoly_init(q2, ctx);
416      fmpq_mpoly_div(res->den, x->den, gd, ctx);
417      fmpq_mpoly_div(q2, y->den, gd, ctx);
418      fmpq_mpoly_mul(res->num, q2, x->num, ctx);
419      fmpq_mpoly_mul(res->den, res->den, y->num, ctx);
420      fmpq_mpoly_add(res->num, res->num, res->den, ctx);
421      fmpq_mpoly_gcd(res->den, res->num, gd, ctx);
422      if (fmpq_mpoly_is_one(res->den, ctx))
423      {
424        fmpq_mpoly_mul(res->den, q2, x->den, ctx);
425      }
426      else
427      {
428        fmpq_mpoly_div(res->num, res->num, res->den, ctx);
429        fmpq_mpoly_div(gd, x->den, res->den, ctx);
430        fmpq_mpoly_mul(res->den, gd, q2, ctx);
431      }
432      fmpq_mpoly_clear(q2, ctx);
433    }
434    fmpq_mpoly_clear(gd, ctx);
435  }
436  #ifdef QA_DEBUG
437  res->p=n_Add(x->p,y->p, ((data_ptr)c->data)->C);
438  #endif
439  n_Test((number)res, c);
440  return (number) res;
441}
442
443static number Div(number a, number b, const coeffs c)
444{
445  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
446  const fmpq_rat_ptr y = (fmpq_rat_ptr) b;
447  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
448  if (fmpq_mpoly_is_zero(y->num, ctx))
449  {
450    WerrorS(nDivBy0);
451    return NULL;
452  }
453  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
454  fmpq_rat_init(res, c);
455  if (fmpq_mpoly_equal(x->den, y->num, ctx)) /* denominators equal */
456  {
457    fmpq_mpoly_mul(res->num, x->num, y->den, ctx);
458    fmpq_mpoly_mul(res->den, x->den, y->num, ctx);
459  }
460  else if (fmpq_mpoly_is_one(x->den, ctx)) /* first denominator 1 */
461  {
462    fmpq_mpoly_t gd;
463    fmpq_mpoly_init(gd, ctx);
464    fmpq_mpoly_gcd(gd, x->num, y->num, ctx);
465    if (fmpq_mpoly_is_one(gd, ctx))
466    {
467      fmpq_mpoly_mul(res->num, x->num, y->den, ctx);
468      fmpq_mpoly_set(res->den, y->num, ctx);
469    }
470    else
471    {
472      fmpq_mpoly_div(res->num, x->num, gd, ctx);
473      fmpq_mpoly_mul(res->num, res->num, y->den, ctx);
474      fmpq_mpoly_div(res->den, y->num, gd, ctx);
475    }
476    fmpq_mpoly_clear(gd, ctx);
477  }
478  else if (fmpq_mpoly_is_one(y->num, ctx)) /* second denominator 1 */
479  {
480    fmpq_mpoly_t gd;
481    fmpq_mpoly_init(gd, ctx);
482    fmpq_mpoly_gcd(gd, y->den, x->den, ctx);
483    if (fmpq_mpoly_is_one(gd, ctx))
484    {
485      fmpq_mpoly_mul(res->num, y->den, x->num, ctx);
486      fmpq_mpoly_set(res->den, x->den, ctx);
487    }
488    else
489    {
490      fmpq_mpoly_div(res->num, y->den, gd, ctx);
491      fmpq_mpoly_mul(res->num, res->num, x->num, ctx);
492      fmpq_mpoly_div(res->den, x->den, gd, ctx);
493    }
494    fmpq_mpoly_clear(gd, ctx);
495  }
496  else /* general case */
497  {
498    fmpq_mpoly_t g1, g2;
499    fmpq_mpoly_ptr n1, n2, d1, d2;
500    fmpq_mpoly_init(g1, ctx);
501    fmpq_mpoly_init(g2, ctx);
502    fmpq_mpoly_gcd(g1, x->num, y->num, ctx);
503    fmpq_mpoly_gcd(g2, y->den, x->den, ctx);
504    n1 = x->num; d2 = y->num;
505    d1 = x->den; n2 = y->den;
506    if (!fmpq_mpoly_is_one(g1, ctx))
507    {
508      fmpq_mpoly_div(res->num, x->num, g1, ctx);
509      fmpq_mpoly_div(g1, y->num, g1, ctx);
510      n1 = res->num; d2 = g1;
511    }
512    if (!fmpq_mpoly_is_one(g2, ctx))
513    {
514      fmpq_mpoly_div(res->den, y->den, g2, ctx);
515      fmpq_mpoly_div(g2, x->den, g2, ctx);
516      n2 = res->den; d1 = g2;
517    }
518    fmpq_mpoly_mul(res->num, n1, n2, ctx);
519    fmpq_mpoly_mul(res->den, d1, d2, ctx);
520    fmpq_mpoly_clear(g1, ctx);
521    fmpq_mpoly_clear(g2, ctx);
522  }
523  fmpq_rat_canonicalise(res, c);
524  #ifdef QA_DEBUG
525  res->p=n_Div(x->p,y->p, ((data_ptr)c->data)->C);
526  #endif
527  n_Test((number)res, c);
528  return (number) res;
529}
530
531static number ExactDiv(number a, number b, const coeffs c)
532{
533  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
534  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
535  const fmpq_rat_ptr y = (fmpq_rat_ptr) b;
536  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
537  if (fmpq_mpoly_is_zero(y->num, ctx))
538  {
539     WerrorS(nDivBy0);
540     return NULL;
541  }
542  fmpq_rat_init(res, c);
543  fmpq_mpoly_div(res->num, x->num, y->num, ctx);
544  assume(fmpq_mpoly_is_one(x->den, ctx));
545  assume(fmpq_mpoly_is_one(y->den, ctx));
546  #ifdef QA_DEBUG
547  res->p=n_ExactDiv(x->p,y->p, ((data_ptr)c->data)->C);
548  #endif
549  n_Test((number)res,c);
550  return (number) res;
551}
552
553// static number IntMod(number a, number b, const coeffs c);
554// {
555// }
556
557static number Init(long i, const coeffs c)
558{
559  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
560  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
561  fmpq_rat_init(res, c);
562  fmpq_mpoly_set_si(res->num, (slong) i, ctx);
563  fmpq_mpoly_set_si(res->den, (slong) 1, ctx);
564  #ifdef QA_DEBUG
565  res->p=n_Init(i, ((data_ptr)c->data)->C);
566  #endif
567  n_Test((number)res,c);
568  return (number) res;
569}
570
571static number InitMPZ(mpz_t i, const coeffs c)
572{
573  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
574  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
575  fmpz_t t;
576  fmpz_init(t);
577  fmpz_set_mpz(t, i);
578  fmpq_rat_init(res, c);
579  fmpq_mpoly_set_fmpz(res->num, t, ctx);
580  fmpq_mpoly_set_si(res->den, (slong) 1, ctx);
581  fmpz_clear(t);
582  #ifdef QA_DEBUG
583  res->p=n_InitMPZ(i, ((data_ptr)c->data)->C);
584  #endif
585  n_Test((number)res,c);
586  return (number) res;
587}
588
589static int Size(number n, const coeffs c)
590{
591  const fmpq_rat_ptr x = (fmpq_rat_ptr) n;
592  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
593  if (fmpq_mpoly_is_zero(x->num, ctx))
594    return 0;
595  unsigned long len=fmpq_mpoly_length(x->num, ctx) +
596         fmpq_mpoly_length(x->den, ctx)-fmpq_mpoly_is_one(x->den, ctx);
597  unsigned long numDegree=fmpq_mpoly_total_degree_si(x->num, ctx);
598  unsigned long denDegree=fmpq_mpoly_total_degree_si(x->den, ctx);
599  unsigned long t= ((numDegree + denDegree)*(numDegree + denDegree) + 1) * len;
600  if (t>INT_MAX) return INT_MAX;
601  else return (int)t;
602}
603
604static long Int(number &n, const coeffs c)
605{
606  const fmpq_rat_ptr x = (fmpq_rat_ptr) n;
607  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
608  if (fmpq_mpoly_is_fmpq(x->den, ctx) && fmpq_mpoly_is_fmpq(x->num, ctx))
609  {
610    long nl = 0;
611    fmpq_t r;
612    fmpq_init(r);
613    fmpq_div(r, x->num->content, x->den->content);
614    if (fmpz_is_one(fmpq_denref(r)))
615    {
616      if (fmpz_fits_si(fmpq_numref(r)))
617        nl = fmpz_get_si(fmpq_numref(r));
618    }
619    fmpq_clear(r);
620    return nl;
621  }
622  return 0;
623}
624
625static void MPZ(mpz_t result, number &n, const coeffs c)
626{
627  mpz_init(result);
628  const fmpq_rat_ptr x = (fmpq_rat_ptr) n;
629  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
630  if (fmpq_mpoly_is_fmpq(x->den, ctx) && fmpq_mpoly_is_fmpq(x->num, ctx))
631  {
632    long nl = 0;
633    fmpq_t r;
634    fmpq_init(r);
635    fmpq_div(r, x->num->content, x->den->content);
636    if (fmpz_is_one(fmpq_denref(r)))
637    {
638      fmpz_get_mpz(result, fmpq_numref(r));
639    }
640    fmpq_clear(r);
641  }
642}
643
644static number Neg(number a, const coeffs c)
645{
646  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
647  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
648  fmpq_mpoly_neg(x->num, x->num, ctx);
649  #ifdef QA_DEBUG
650  x->p=n_InpNeg(x->p, ((data_ptr)c->data)->C);
651  #endif
652  n_Test((number)a, c);
653  return a;
654}
655
656static number Invers(number a, const coeffs c)
657{
658  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
659  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
660  if (fmpq_mpoly_is_zero(x->num, ctx))
661  {
662    WerrorS(nDivBy0);
663    return NULL;
664  }
665  else
666  {
667    fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
668    fmpq_rat_init(res, c);
669    fmpq_mpoly_set(res->num, x->den, ctx);
670    fmpq_mpoly_set(res->den, x->num, ctx);
671  #ifdef QA_DEBUG
672  res->p=n_Invers(x->p, ((data_ptr)c->data)->C);
673  #endif
674  n_Test((number)res, c);
675    return (number) res;
676  }
677}
678
679static number Copy(number a, const coeffs c)
680{
681  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
682  fmpq_rat_init(res, c);
683  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
684  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
685  fmpq_mpoly_set(res->num, x->num, ctx);
686  fmpq_mpoly_set(res->den, x->den, ctx);
687  #ifdef QA_DEBUG
688  res->p=n_Copy(x->p, ((data_ptr)c->data)->C);
689  #endif
690  n_Test((number)res, c);
691  return (number) res;
692}
693
694//static number RePart(number a, const coeffs c)
695//{
696//}
697
698//static number ImPart(number a, const coeffs c)
699//{
700//}
701
702static BOOLEAN IsOne(number a, const coeffs c)
703{
704  if (a==NULL) return FALSE;
705  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
706  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
707  if (!fmpq_mpoly_is_fmpq(x->num, ctx))
708    return FALSE;
709  if (!fmpq_mpoly_is_fmpq(x->den, ctx))
710    return FALSE;
711  return fmpq_equal(x->num->content, x->den->content);
712}
713
714static BOOLEAN IsZero(number a, const coeffs c)
715{
716  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
717  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
718  return fmpq_mpoly_is_zero(x->num, ctx);
719}
720
721static void WriteLong(number a, const coeffs c)
722{
723  if (a==NULL)
724  {
725    StringAppendS("o");
726    return;
727  }
728  n_Test(a,c);
729  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
730  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
731  if (fmpq_mpoly_is_zero(x->den, ctx))
732  {
733    StringAppendS("?/o");
734    return;
735  }
736  fmpz_t t;
737  long int i, j, k, nmax_i, dmax_i, max_digits;
738  fmpq_rat_canonicalise(x, c);
739  if (fmpq_mpoly_is_zero(x->num, ctx))
740    StringAppendS("0");
741  else
742  {
743    BOOLEAN num_is_const = fmpq_mpoly_is_fmpq(x->num, ctx);
744    BOOLEAN den_is_const = fmpq_mpoly_is_fmpq(x->den, ctx);
745    BOOLEAN need_times;
746    fmpq_mpoly_struct * znum = x->num;
747    fmpq_mpoly_struct * zden = x->den;
748    slong nvars = fmpq_mpoly_ctx_nvars(ctx);
749    fmpz_init(t);
750    nmax_i = 0;
751    dmax_i = 0;
752    for (i = 1; i < fmpq_mpoly_length(znum, ctx); i++)
753    {
754      if (fmpz_cmpabs(fmpq_mpoly_zpoly_term_coeff_ref(znum, i, ctx),
755                      fmpq_mpoly_zpoly_term_coeff_ref(znum, nmax_i, ctx)) > 0)
756      {
757         nmax_i = i;
758      }
759    }
760    for (i = 1; i < fmpq_mpoly_length(zden, ctx); i++)
761    {
762      if (fmpz_cmpabs(fmpq_mpoly_zpoly_term_coeff_ref(zden, i, ctx),
763                      fmpq_mpoly_zpoly_term_coeff_ref(zden, dmax_i, ctx)) > 0)
764      {
765        dmax_i = i;
766      }
767    }
768    if (fmpz_cmpabs(fmpq_mpoly_zpoly_term_coeff_ref(znum, nmax_i, ctx),
769                       fmpq_mpoly_zpoly_term_coeff_ref(zden, dmax_i, ctx)) > 0)
770    {
771      fmpz_mul(t, fmpq_numref(x->num->content),
772                  fmpq_mpoly_zpoly_term_coeff_ref(znum, nmax_i, ctx));
773      max_digits = fmpz_sizeinbase(t, 10);
774    }
775    else
776    {
777      fmpz_mul(t, fmpq_numref(x->den->content),
778                  fmpq_mpoly_zpoly_term_coeff_ref(zden, dmax_i, ctx));
779      max_digits =fmpz_sizeinbase(t, 10);
780    }
781    char *s = (char*) omAlloc(max_digits + 5);
782    if (!num_is_const)
783      StringAppendS("(");
784    if (fmpq_mpoly_is_one(x->num, ctx))
785      StringAppendS("1");
786    else
787    {
788      for (i = 0; i < fmpq_mpoly_length(x->num, ctx); i++)
789      {
790        need_times = TRUE;
791        fmpz_mul(t, fmpq_mpoly_zpoly_term_coeff_ref(znum, i, ctx),
792                                       fmpq_numref(x->num->content));
793        if (i != 0 && fmpz_sgn(t) > 0)
794          StringAppendS("+");
795        BOOLEAN need_1=FALSE;
796        if (!fmpz_is_one(t))
797        {
798          fmpz_get_str(s, 10, t);
799          {
800            int l=strlen(s);
801            while((l>0)&&(!isdigit(s[l]))) l--;
802            s[l+1]='\0';
803          }
804          if (strcmp(s,"-1")==0)
805          {
806            StringAppendS("-");
807            need_1 = TRUE;
808            need_times = FALSE;
809          }
810          else
811            StringAppendS(s);
812        }
813        else
814        {
815          need_1 = TRUE;
816          need_times = FALSE;
817        }
818        for (j = 0; j < c->iNumberOfParameters; j++)
819        {
820          k = fmpq_mpoly_get_term_var_exp_ui(x->num, i, j, ctx);
821          if (k != 0)
822          {
823            need_1 = FALSE;
824            if (need_times)
825              StringAppendS("*");
826            if (k != 1)
827              StringAppend("%s^%d", c->pParameterNames[j], k);
828            else
829              StringAppendS(c->pParameterNames[j]);
830            need_times = TRUE;
831          }
832        }
833        if (need_1)  StringAppendS("1");
834      }
835    }
836    if (!num_is_const)
837       StringAppendS(")");
838    if (!fmpq_mpoly_is_one(x->den, ctx))
839    {
840      BOOLEAN closing_paren=FALSE;
841      StringAppendS("/");
842      if (!den_is_const)
843      {
844        StringAppendS("(");
845        closing_paren = TRUE;
846      }
847      for (i = 0; i < fmpq_mpoly_length(x->den, ctx); i++)
848      {
849        need_times = TRUE;
850        fmpz_mul(t, fmpq_mpoly_zpoly_term_coeff_ref(zden, i, ctx),
851                                       fmpq_numref(x->den->content));
852        if (i == 0)
853        {
854          if ((fmpz_sgn(t) < 0) && den_is_const)
855          {
856            StringAppendS("(");
857            closing_paren = TRUE;
858          }
859        }
860        else if (fmpz_sgn(t) > 0)
861          StringAppendS("+");
862        if (!fmpz_is_one(t))
863        {
864          fmpz_get_str(s, 10, t);
865          {
866            int l=strlen(s);
867            while((l>0)&&(!isdigit(s[l]))) l--;
868            s[l+1]='\0';
869          }
870          StringAppendS(s);
871        }
872        else
873        {
874          need_times = FALSE;
875        }
876        for (j = 0; j < nvars; j++)
877        {
878          k = fmpq_mpoly_get_term_var_exp_ui(x->den, i, j, ctx);
879          if (k != 0)
880          {
881            if (need_times)
882              StringAppendS("*");
883            if (k != 1)
884              StringAppend("%s^%d", c->pParameterNames[j], k);
885            else
886              StringAppendS(c->pParameterNames[j]);
887            need_times = TRUE;
888          }
889        }
890      }
891      if (closing_paren)
892        StringAppendS(")");
893    }
894    fmpz_clear(t);
895    omFree(s);
896  }
897}
898
899static void WriteShort(number a, const coeffs c)
900{
901  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
902  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
903  fmpz_t t;
904  char *s;
905  long int i, j, k, nmax_i, dmax_i, max_digits;
906  fmpq_rat_canonicalise(x, c);
907  if (fmpq_mpoly_is_zero(x->num, ctx))
908    StringAppendS("0");
909  else
910  {
911    BOOLEAN num_is_const = fmpq_mpoly_is_fmpq(x->num, ctx);
912    BOOLEAN den_is_const = fmpq_mpoly_is_fmpq(x->den, ctx);
913    fmpq_mpoly_struct * znum = x->num;
914    fmpq_mpoly_struct * zden = x->den;
915    slong nvars = fmpq_mpoly_ctx_nvars(ctx);
916    fmpz_init(t);
917    nmax_i = 0;
918    dmax_i = 0;
919    for (i = 1; i < fmpq_mpoly_length(znum, ctx); i++)
920    {
921      if (fmpz_cmpabs(fmpq_mpoly_zpoly_term_coeff_ref(znum, i, ctx),
922                      fmpq_mpoly_zpoly_term_coeff_ref(znum, nmax_i, ctx)) > 0)
923      {
924         nmax_i = i;
925      }
926    }
927    for (i = 1; i < fmpq_mpoly_length(zden, ctx); i++)
928    {
929      if (fmpz_cmpabs(fmpq_mpoly_zpoly_term_coeff_ref(zden, i, ctx),
930                      fmpq_mpoly_zpoly_term_coeff_ref(zden, dmax_i, ctx)) > 0)
931      {
932        dmax_i = i;
933      }
934    }
935    if (fmpz_cmpabs(fmpq_mpoly_zpoly_term_coeff_ref(znum, nmax_i, ctx),
936                       fmpq_mpoly_zpoly_term_coeff_ref(zden, dmax_i, ctx)) > 0)
937    {
938      fmpz_mul(t, fmpq_numref(x->num->content),
939                  fmpq_mpoly_zpoly_term_coeff_ref(znum, nmax_i, ctx));
940      max_digits = fmpz_sizeinbase(t, 10);
941    } else
942    {
943      fmpz_mul(t, fmpq_numref(x->den->content),
944                  fmpq_mpoly_zpoly_term_coeff_ref(zden, dmax_i, ctx));
945      max_digits = fmpz_sizeinbase(t, 10);
946    }
947    s = (char*) omAlloc(max_digits + 2);
948    if (!num_is_const)
949      StringAppendS("(");
950    if (fmpq_mpoly_is_one(x->num, ctx))
951      StringAppendS("1");
952    else
953    {
954      for (i = 0; i < fmpq_mpoly_length(x->num, ctx); i++)
955      {
956        fmpz_mul(t, fmpq_mpoly_zpoly_term_coeff_ref(znum, i, ctx),
957                                       fmpq_numref(x->num->content));
958        if (i != 0 && fmpz_sgn(t) > 0)
959          StringAppendS("+");
960        if (!fmpz_is_one(t))
961        {
962          fmpz_get_str(s, 10, t);
963          StringAppendS(s);
964        }
965        for (j = 0; j < nvars; j++)
966        {
967          k = fmpq_mpoly_get_term_var_exp_ui(x->num, i, j, ctx);
968          if (k != 0)
969          {
970            if (k != 1)
971              StringAppend("%s%d", c->pParameterNames[j], k);
972            else
973              StringAppendS(c->pParameterNames[j]);
974          }
975        }
976      }
977    }
978    if (!num_is_const)
979       StringAppendS(")");
980    if (!fmpq_mpoly_is_one(x->den, ctx))
981    {
982      StringAppendS("/");
983      if (!den_is_const)
984        StringAppendS("(");
985      for (i = 0; i < fmpq_mpoly_length(x->den, ctx); i++)
986      {
987        fmpz_mul(t, fmpq_mpoly_zpoly_term_coeff_ref(zden, i, ctx),
988                                       fmpq_numref(x->den->content));
989        if (i != 0 && fmpz_sgn(t) > 0)
990          StringAppendS("+");
991        if (!fmpz_is_one(t))
992        {
993          fmpz_get_str(s, 10, t);
994          StringAppendS(s);
995        }
996        for (j = 0; j < nvars; j++)
997        {
998          k = fmpq_mpoly_get_term_var_exp_ui(x->num, i, j, ctx);
999          if (k != 0)
1000          {
1001            if (k != 1)
1002              StringAppend("%s%d", c->pParameterNames[j], k);
1003            else
1004              StringAppendS(c->pParameterNames[j]);
1005          }
1006        }
1007      }
1008      if (!den_is_const)
1009        StringAppendS(")");
1010    }
1011    fmpz_clear(t);
1012  }
1013}
1014
1015static const char* Read(const char * st, number * a, const coeffs c)
1016{
1017  // we only read "monomials" (i.e. [-][digits][parameter]),
1018  // everything else (+,*,^,()) is left to the singular interpreter
1019  long int j;
1020  char *s = (char *) st;
1021  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1022  slong nvars = fmpq_mpoly_ctx_nvars(ctx);
1023  *a = (number) omAlloc(sizeof(fmpq_rat_struct));
1024  fmpq_rat_init((fmpq_rat_ptr)(*a), c);
1025  BOOLEAN neg = FALSE;
1026  if (*s=='-')
1027  {
1028    neg = TRUE;
1029    s++;
1030  }
1031  if (isdigit(*s))
1032  {
1033    fmpz_t z;
1034    fmpz_init(z);
1035    s = nlEatLong((char *) s, z);
1036    fmpq_mpoly_set_fmpz(((fmpq_rat_ptr)(*a))->num, z, ctx);
1037    fmpq_mpoly_one(((fmpq_rat_ptr)(*a))->den, ctx);
1038    if (*s == '/')
1039    {
1040      s++;
1041      s = nlEatLong((char *) s, z);
1042      fmpq_mpoly_scalar_div_fmpz(((fmpq_rat_ptr)(*a))->num,
1043                                    ((fmpq_rat_ptr)(*a))->num, z, ctx);
1044    }
1045    fmpz_clear(z);
1046  }
1047  else
1048  {
1049    BOOLEAN found=FALSE;
1050    for (j = 0; j < nvars; j++)
1051    {
1052      if (strncmp(s, c->pParameterNames[j],
1053                                  strlen(c->pParameterNames[j])) == 0)
1054      {
1055        found=TRUE;
1056        fmpq_mpoly_gen(((fmpq_rat_ptr)(*a))->num, j, ctx);
1057        s += strlen(c->pParameterNames[j]);
1058        if (isdigit(*s))
1059        {
1060          int i = 1;
1061          s = nEati(s, &i, 0);
1062          if (i != 1)
1063          {
1064            fmpq_mpoly_pow_ui(((fmpq_rat_ptr)(*a))->num,
1065                         ((fmpq_rat_ptr)(*a))->num, (long int) i, ctx);
1066          }
1067        }
1068      }
1069    }
1070    if (!found) fmpq_mpoly_one(((fmpq_rat_ptr)(*a))->num, ctx);
1071    fmpq_mpoly_one(((fmpq_rat_ptr)(*a))->den, ctx);
1072  }
1073  if (neg)
1074    fmpq_mpoly_neg(((fmpq_rat_ptr)(*a))->num, ((fmpq_rat_ptr)(*a))->num, ctx);
1075  #ifdef QA_DEBUG
1076  poly pp=convFlintMPSingP(((fmpq_rat_ptr)(*a))->num,ctx,((data_ptr)c->data)->C->extRing);
1077  fraction f=(fraction)n_Init(1,((data_ptr)c->data)->C); /*leak*/
1078  NUM(f)=pp;
1079  ((fmpq_rat_ptr)(*a))->p=(number)f;
1080  #endif
1081  n_Test((*a),c);
1082  return s;
1083}
1084
1085static BOOLEAN Greater(number a, number b, const coeffs c)
1086{
1087  return Size(a, c) > Size(b, c);
1088}
1089
1090static void Delete(number * a, const coeffs c)
1091{
1092  if ((*a) != NULL)
1093  {
1094    const fmpq_rat_ptr x = (fmpq_rat_ptr) *a;
1095    fmpq_rat_clear(x, c);
1096    #ifdef QA_DEBUG
1097    n_Delete(&(x->p),((data_ptr)c->data)->C);
1098    #endif
1099    omFree(*a);
1100    *a = NULL;
1101  }
1102}
1103
1104static BOOLEAN Equal(number a, number b, const coeffs c)
1105{
1106  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
1107  const fmpq_rat_ptr y = (fmpq_rat_ptr) b;
1108  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1109  if (!fmpz_mpoly_equal(x->num->zpoly, y->num->zpoly, ctx->zctx))
1110  {
1111    return FALSE;
1112  }
1113  if (!fmpz_mpoly_equal(x->den->zpoly, y->den->zpoly, ctx->zctx))
1114  {
1115    return FALSE;
1116  }
1117  fmpz_t t1, t2;
1118  fmpz_init(t1);
1119  fmpz_init(t2);
1120  fmpz_mul(t1, fmpq_numref(x->num->content), fmpq_denref(x->den->content));
1121  fmpz_mul(t1, t1, fmpq_denref(y->num->content));
1122  fmpz_mul(t1, t1, fmpq_numref(y->den->content));
1123  fmpz_mul(t2, fmpq_numref(y->num->content), fmpq_denref(y->den->content));
1124  fmpz_mul(t2, t2, fmpq_denref(x->num->content));
1125  fmpz_mul(t2, t2, fmpq_numref(x->den->content));
1126  int eq = fmpz_equal(t1, t2);
1127  fmpz_clear(t1);
1128  fmpz_clear(t2);
1129  return eq;
1130}
1131
1132static BOOLEAN IsMOne(number a, const coeffs c)
1133{
1134  if (a==NULL) return FALSE;
1135  fmpq_t content;
1136  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
1137  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1138  if (!fmpq_mpoly_is_fmpq(x->num, ctx))
1139    return FALSE;
1140  if (!fmpq_mpoly_is_fmpq(x->den, ctx))
1141    return FALSE;
1142  fmpq_init(content);
1143  fmpq_neg(content, x->num->content);
1144  int eq = fmpq_equal(content, x->den->content);
1145  fmpq_clear(content);
1146  return eq;
1147}
1148
1149static BOOLEAN GreaterZero(number, const coeffs)
1150{
1151  return TRUE; /* everything in parens for now so need + sign */
1152}
1153
1154static void Power(number a, int i, number * result, const coeffs c)
1155{
1156  *result= (number) omAlloc(sizeof(fmpq_rat_struct));
1157  fmpq_rat_init((fmpq_rat_ptr) (*result), c);
1158  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
1159  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1160  fmpq_mpoly_pow_ui(((fmpq_rat_ptr)(*result))->num, x->num, (slong) i, ctx);
1161  fmpq_mpoly_pow_ui(((fmpq_rat_ptr)(*result))->den, x->den, (slong) i, ctx);
1162}
1163
1164static number GetDenom(number &n, const coeffs c)
1165{
1166  const fmpq_rat_ptr x = (fmpq_rat_ptr) n;
1167  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1168  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
1169  fmpq_rat_init(res, c);
1170  fmpq_mpoly_set(res->num, x->den, ctx);
1171  fmpq_mpoly_one(res->den, ctx);
1172  return (number) res;
1173}
1174
1175static number GetNumerator(number &n, const coeffs c)
1176{
1177  const fmpq_rat_ptr x = (fmpq_rat_ptr) n;
1178  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1179  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
1180  fmpq_rat_init(res, c);
1181  fmpq_mpoly_set(res->num, x->num, ctx);
1182  fmpq_mpoly_one(res->den, ctx);
1183  return (number) res;
1184}
1185
1186static number ExtGcd(number a, number b, number *s, number *t, const coeffs c)
1187{
1188  WerrorS("not a Euclidean ring: ExtGcd");
1189  return NULL;
1190}
1191
1192static number Lcm(number a, number b, const coeffs c)
1193{
1194  WerrorS("not yet: Lcm");
1195  return NULL;
1196}
1197
1198static number Q2Frac(number a, const coeffs src, const coeffs dst)
1199{
1200  number res;
1201  if (SR_HDL(a) & SR_INT)
1202  {
1203     res=Init(SR_TO_INT(a),dst);
1204     n_Test(res,dst);
1205  }
1206  else if (a->s==3)
1207  {
1208    res=InitMPZ(a->z,dst);
1209    n_Test(res,dst);
1210  }
1211  else
1212  {
1213    number z=InitMPZ(a->z,dst);
1214    number n=InitMPZ(a->n,dst);
1215    res=Div(z,n,dst);
1216    Delete(&z,dst);
1217    Delete(&n,dst);
1218    n_Test(res,dst);
1219    return res;
1220  }
1221  return res;
1222}
1223
1224static number Z2Frac(number a, const coeffs src, const coeffs dst)
1225{
1226  return InitMPZ((mpz_ptr)a,dst);
1227}
1228
1229static number Zp2Frac(number a, const coeffs src, const coeffs dst)
1230{
1231  return Init(n_Int(a,src),dst);
1232}
1233
1234static nMapFunc SetMap(const coeffs src, const coeffs dst)
1235{
1236  if (src == dst) return ndCopyMap;
1237  if (nCoeff_is_Q_or_BI(src) && (src->rep==n_rep_gap_rat))  /*Q, coeffs_BIGINT */
1238    return Q2Frac;
1239  if(src->rep==n_rep_gap_gmp) /*Z */
1240    return Z2Frac;
1241  if(nCoeff_is_Zp(src))
1242    return Zp2Frac;
1243
1244  return NULL;
1245}
1246
1247//static void InpMult(number &a, number b, const coeffs c)
1248//{
1249//}
1250
1251//static void InpAdd(number &a, number b, const coeffs c)
1252//{
1253//}
1254
1255#if 0
1256static number Init_bigint(number i, const coeffs dummy, const coeffs dst)
1257{
1258  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)(dst->data))->ctx;
1259  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
1260  fmpz_t f;
1261  fmpq_rat_init(res, dst);
1262  if (SR_HDL(i) & SR_INT)
1263  {
1264    fmpq_mpoly_set_si(res->num, SR_TO_INT(i), ctx);
1265  }
1266  else
1267  {
1268    fmpz_init(f);
1269    fmpz_set_mpz(f, i->z);
1270    fmpq_mpoly_set_fmpz(res->num, f, ctx);
1271    fmpz_clear(f);
1272  }
1273  fmpq_mpoly_set_si(res->den, 1, ctx);
1274  return (number) res;
1275}
1276#endif
1277
1278#if 0
1279static number Farey(number p, number n, const coeffs c)
1280{
1281  WerrorS("not yet: Farey");
1282  return NULL;
1283}
1284#endif
1285
1286#if 0
1287static number ChineseRemainder(number *x, number *q, int rl,
1288                    BOOLEAN sym, CFArray &inv_cache, const coeffs c)
1289{
1290  WerrorS("not yet: ChineseRemainder");
1291  return NULL;
1292}
1293#endif
1294
1295static int ParDeg(number a, const coeffs c)
1296{
1297    const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
1298    const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1299    return (int) (fmpq_mpoly_total_degree_si(x->num, ctx) -
1300                  fmpq_mpoly_total_degree_si(x->den, ctx));
1301}
1302
1303static number Parameter(const int i, const coeffs c)
1304{
1305  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1306  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
1307  fmpq_rat_init(res, c);
1308  fmpq_mpoly_gen(res->num, (slong) i, ctx);
1309  fmpq_mpoly_one(res->den, ctx);
1310  return (number) res;
1311}
1312
1313static number SubringGcd(number a, number b, const coeffs c)
1314{
1315  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
1316  fmpq_rat_init(res, c);
1317  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
1318  const fmpq_rat_ptr y = (fmpq_rat_ptr) b;
1319  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1320  fmpq_mpoly_gcd(res->num, x->num, y->num, ctx);
1321  // handle content:
1322  fmpz_t cont;
1323  fmpz_init(cont);
1324  fmpz_gcd(cont, fmpq_numref(x->num->content), fmpq_numref(y->num->content));
1325  if (!fmpz_is_one(cont))
1326  {
1327     fmpq_mul_fmpz(res->num->content, res->num->content, cont);
1328  }
1329  fmpz_gcd(cont, fmpq_denref(x->num->content), fmpq_denref(y->num->content));
1330  if (!fmpz_is_one(cont))
1331  {
1332     fmpq_div_fmpz(res->num->content, res->num->content, cont);
1333  }
1334  fmpz_clear(cont);
1335  fmpq_mpoly_one(res->den, ctx);
1336  fmpq_rat_canonicalise(res, c);
1337  #ifdef QA_DEBUG
1338  res->p=n_SubringGcd(x->p,y->p, ((data_ptr)c->data)->C);
1339  #endif
1340  n_Test((number)res, c);
1341  return (number) res;
1342}
1343
1344static number NormalizeHelper(number a, number b, const coeffs c)
1345{
1346  fmpq_rat_ptr res = (fmpq_rat_ptr) omAlloc(sizeof(fmpq_rat_struct));
1347  fmpq_rat_init(res, c);
1348  const fmpq_rat_ptr x = (fmpq_rat_ptr) a;
1349  const fmpq_rat_ptr y = (fmpq_rat_ptr) b;
1350  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1351  fmpq_mpoly_t gd;
1352  fmpq_mpoly_init(gd,ctx);
1353  fmpq_mpoly_one(gd,ctx); // value for gd, if fmpq_mpoly_gcd fails
1354  fmpq_mpoly_gcd(gd, x->num, y->den, ctx);
1355  fmpq_mpoly_mul(res->num,  x->num, y->den, ctx);
1356  if (!fmpq_mpoly_is_one(gd, ctx))// &&(!fmpq_mpoly_is_zero(gd, ctx)))
1357    fmpq_mpoly_div(res->num, res->num, gd, ctx);
1358  fmpq_mpoly_one(res->den, ctx);
1359  #ifdef QA_DEBUG
1360  res->p=n_NormalizeHelper(x->p,y->p, ((data_ptr)c->data)->C);
1361  #endif
1362  n_Test((number)res, c);
1363  return (number) res;
1364}
1365
1366#if 0
1367static void WriteFd(number a, const ssiInfo *d, const coeffs c)
1368{
1369  // format: len a_len(num den) .. a_0
1370/*  Currently not implemented
1371  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)c->data)->ctx;
1372  const fmpq_rat_ptr aa = (fmpq_rat_ptr) a;
1373  int l = fmpq_mpoly_length(aa->num, ctx);
1374  fprintf(d->f_write, "%d ", l);
1375  mpq_t m;
1376  mpq_init(m);
1377  mpz_t num, den;
1378  mpz_init(num);
1379  mpz_init(den);
1380  for(int i = l; i >= 0; i--)
1381  {
1382    fmpq_mpoly_get_coeff_mpq(m, aa->num, i);
1383    mpq_get_num(num, m);
1384    mpq_get_den(den, m);
1385    mpz_out_str(d->f_write, SSI_BASE, num);
1386    fputc(' ', d->f_write);
1387    mpz_out_str(d->f_write, SSI_BASE, den);
1388    fputc(' ', d->f_write);
1389  }
1390  mpz_clear(den);
1391  mpz_clear(num);
1392  mpq_clear(m);
1393*/
1394}
1395#endif
1396
1397#if 0
1398static number ReadFd(const ssiInfo *d, const coeffs c)
1399{
1400  // format: len a_len .. a_0
1401/* Currently not implemented
1402  fmpq_mpoly_ptr aa = (fmpq_mpoly_ptr) omAlloc(sizeof(fmpq_mpoly_t));
1403  fmpq_mpoly_init(aa);
1404  int l = s_readint(d->f_read);
1405  mpz_t nm;
1406  mpz_init(nm);
1407  mpq_t m;
1408  mpq_init(m);
1409  for (int i = l; i >= 0; i--)
1410  {
1411    s_readmpz_base(d->f_read, nm, SSI_BASE);
1412    mpq_set_num(m, nm);
1413    s_readmpz_base(d->f_read, nm, SSI_BASE);
1414    mpq_set_den(m, nm);
1415    fmpq_mpoly_set_coeff_mpq(aa, i, m);
1416  }
1417  mpz_clear(nm);
1418  mpq_clear(m);
1419  return (number)aa;
1420*/
1421  return NULL;
1422}
1423#endif
1424
1425// cfClearContent
1426
1427// cfClearDenominators
1428
1429#if 0
1430static number ConvFactoryNSingN(const CanonicalForm n, const coeffs c)
1431{
1432  WerrorS("not yet: ConvFactoryNSingN");
1433  return NULL;
1434}
1435#endif
1436
1437#if 0
1438static CanonicalForm ConvSingNFactoryN(number n, BOOLEAN setChar, const coeffs c)
1439{
1440  WerrorS("not yet: ConvSingNFactoryN");
1441  return CanonicalForm(0);
1442}
1443#endif
1444
1445char * QratCoeffName(const coeffs c)
1446{
1447  STATIC_VAR char CoeffName_flint_Qrat[200];
1448  sprintf(CoeffName_flint_Qrat, "flintQQ(%s",c->pParameterNames[0]);
1449  for(int i=1; i<c->iNumberOfParameters;i++)
1450  {
1451    strcat(CoeffName_flint_Qrat,",");
1452    strcat(CoeffName_flint_Qrat,c->pParameterNames[i]);
1453  }
1454  strcat(CoeffName_flint_Qrat,")");
1455  return (char*) CoeffName_flint_Qrat;
1456
1457}
1458
1459coeffs flintQratInitCfByName(char *s, n_coeffType n)
1460{
1461  const char start[] = "flintQ(";
1462  const int start_len = strlen(start);
1463  if (strncmp(s, start, start_len) == 0)
1464  {
1465    s += start_len;
1466    // count ,
1467    char *p=s;
1468    int N=0;
1469    loop
1470    {
1471      while((*p!=',')&&(*p!=')')&&(*p!='\0')) p++;
1472      if (*p==',')       { p++; N++;}
1473      else if (*p==')')  { p++; N++; break;}
1474      else if (*p=='\0') { break;}
1475    }
1476    // get names
1477    char *names[N];
1478    int i=0;
1479    p=s;
1480    loop
1481    {
1482      while((*p!=',')&&(*p!=')')&&(*p!='\0')) p++;
1483      if ((*p==',')||(*p=')'))
1484      {
1485        char c=*p;
1486        *p='\0';
1487        names[i]=omStrDup(s);
1488        *p=c;
1489        i++;
1490        p++;
1491        s=p;
1492        if (c==')') break;
1493      }
1494      if (*p=='\0') break;
1495    }
1496    QaInfo pp;
1497    pp.N=N;
1498    pp.names=names;
1499    coeffs cf=nInitChar(n,&pp);
1500    for(i=0;i<N;i++) omFree(names[i]);
1501    return cf;
1502  }
1503  return NULL;
1504}
1505
1506#ifdef LDEBUG
1507static BOOLEAN DBTest(number c, const char *, const int, const coeffs cf)
1508{
1509  const fmpq_rat_ptr x = (fmpq_rat_ptr) c;
1510  if ((x!=NULL)
1511  && ((x->num==NULL)||(x->den==NULL)))
1512  {
1513    dReportError("NULL num., or den.\n");
1514    return FALSE;
1515  }
1516  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)cf->data)->ctx;
1517  fmpq_mpoly_assert_canonical(x->num,ctx);
1518  fmpq_mpoly_assert_canonical(x->den,ctx);
1519  if (fmpq_mpoly_is_zero(x->den, ctx))
1520  {
1521    dReportError("den.==0\n");
1522    return FALSE;
1523  }
1524  fmpz_t n, d;
1525  fmpz_init(n);
1526  fmpz_init(d);
1527  fmpz_gcd(n, fmpq_numref(x->num->content), fmpq_numref(x->den->content));
1528  fmpz_lcm(d, fmpq_denref(x->num->content), fmpq_denref(x->den->content));
1529  if (!fmpz_is_one(d))
1530  {
1531     dReportError("canon needed (1)");
1532     return TRUE;
1533  }
1534  if (!fmpz_is_one(n))
1535  {
1536     dReportError("canon needed (2)");
1537     return TRUE;
1538  }
1539  fmpz_clear(n);
1540  fmpz_clear(d);
1541  #ifdef QA_DEBUG
1542  poly pp=convFlintMPSingP(x->num,ctx,((data_ptr)cf->data)->C->extRing);
1543  fraction f=(fraction)x->p;
1544  if (f==NULL)
1545  {
1546    dReportError("x->p==NULL\n");
1547    return FALSE;
1548  }
1549  else
1550  {
1551    if (!p_EqualPolys(pp,NUM(f),((data_ptr)cf->data)->C->extRing))
1552    {
1553      p_Write(pp,((data_ptr)cf->data)->C->extRing);
1554      PrintS("num, p=");
1555      p_Write(NUM(f),((data_ptr)cf->data)->C->extRing);
1556      dReportError("num wrong.\n");
1557      return FALSE;
1558    }
1559    if (DEN(f)!=NULL)
1560    {
1561      pp=convFlintMPSingP(x->den,ctx,((data_ptr)cf->data)->C->extRing);
1562      if (!p_EqualPolys(pp,DEN(f),((data_ptr)cf->data)->C->extRing))
1563      {
1564        p_Write(pp,((data_ptr)cf->data)->C->extRing);
1565        PrintS("den, p=");
1566        p_Write(NUM(f),((data_ptr)cf->data)->C->extRing);
1567        dReportError("den wrong.\n");
1568        return FALSE;
1569      }
1570    }
1571  }
1572  #endif
1573  return TRUE;
1574}
1575
1576#endif
1577static void KillChar(coeffs cf)
1578{
1579  for(int i=0;i<cf->iNumberOfParameters;i++)
1580    omFree((ADDRESS)(cf->pParameterNames[i]));
1581  omFreeSize(cf->pParameterNames,sizeof(char*));
1582  const fmpq_ctx_ptr ctx = (fmpq_ctx_ptr) ((data_ptr)cf->data)->ctx;
1583  fmpq_mpoly_ctx_clear(ctx);
1584  omFree(cf->data);
1585}
1586
1587BOOLEAN flintQrat_InitChar(coeffs cf, void * infoStruct)
1588{
1589  QaInfo *pp=(QaInfo*)infoStruct;
1590  cf->cfCoeffName    = QratCoeffName;
1591  cf->nCoeffIsEqual  = CoeffIsEqual;
1592  cf->cfKillChar     = KillChar;
1593  cf->ch             = 0; //char 0
1594  cf->cfMult         = Mult;
1595  cf->cfSub          = Sub;
1596  cf->cfAdd          = Add;
1597  cf->cfDiv          = Div;
1598  cf->cfExactDiv     = Div; // ???
1599  cf->cfInit         = Init;
1600  cf->cfInitMPZ      = InitMPZ;
1601  cf->cfSize         = Size;
1602  cf->cfInt          = Int;
1603  cf->cfMPZ          = MPZ;
1604  cf->cfInpNeg       = Neg;
1605  cf->cfInvers       = Invers;
1606  cf->cfCopy         = Copy;
1607  cf->cfRePart       = Copy;
1608  // default: cf->cfImPart       = ndReturn0;
1609  cf->cfWriteLong    = WriteLong;
1610  cf->cfWriteShort   = WriteLong;
1611  cf->cfRead         = Read;
1612  //cf->cfNormalize    = Normalize;
1613
1614  //cf->cfDivComp=
1615  //cf->cfIsUnit=
1616  //cf->cfGetUnit=
1617  //cf->cfDivBy=
1618
1619  cf->cfGreater      = Greater;
1620  cf->cfEqual        = Equal;
1621  cf->cfIsZero       = IsZero;
1622  cf->cfIsOne        = IsOne;
1623  cf->cfIsMOne       = IsMOne;
1624  cf->cfGreaterZero  = GreaterZero;
1625
1626  cf->cfPower        = Power;
1627  cf->cfGetDenom     = GetDenom;
1628  cf->cfGetNumerator = GetNumerator;
1629  cf->cfExtGcd       = ExtGcd;
1630  cf->cfSubringGcd   = SubringGcd;
1631  cf->cfNormalizeHelper= NormalizeHelper;
1632  cf->cfLcm          = Lcm;
1633  cf->cfDelete       = Delete;
1634  cf->cfSetMap       = SetMap;
1635  // default: cf->cfInpMult
1636  // default: cf->cfInpAdd
1637  //cf->cfFarey        =Farey;
1638  //cf->cfChineseRemainder = ChineseRemainder;
1639  cf->cfParDeg       = ParDeg;
1640  cf->cfParameter    = Parameter;
1641  //  cf->cfClearContent = ClearContent;
1642  //  cf->cfClearDenominators = ClearDenominators;
1643  //cf->convFactoryNSingN = ConvFactoryNSingN;
1644  //cf->convSingNFactoryN = ConvSingNFactoryN;
1645  //cf->cfWriteFd      = WriteFd;
1646  //cf->cfReadFd       = ReadFd;
1647#ifdef LDEBUG
1648  cf->cfDBTest       = DBTest;
1649#endif
1650
1651  cf->iNumberOfParameters = pp->N;
1652  char **pn = (char**) omAlloc0(pp->N*sizeof(char*));
1653  for(int i=0;i<pp->N;i++)
1654  {
1655    pn[i] = omStrDup(pp->names[i]);
1656  }
1657  cf->pParameterNames = (const char **) pn;
1658  cf->has_simple_Inverse = FALSE;
1659  cf->has_simple_Alloc = FALSE;
1660  cf->is_field = TRUE;
1661  cf->is_domain = TRUE;
1662
1663  fmpq_rat_data_struct *ps=(fmpq_rat_data_struct*)omAlloc(sizeof(fmpq_rat_data_struct));
1664  ps->ctx=(fmpq_mpoly_ctx_struct*)omAlloc(sizeof(fmpq_mpoly_ctx_struct));
1665  #ifdef QA_DEBUG
1666  ps->C=pp->C;
1667  #endif
1668  fmpq_mpoly_ctx_init(ps->ctx,pp->N,ORD_LEX);
1669  cf->data=ps;
1670  return FALSE;
1671}
1672#else
1673BOOLEAN flintQrat_InitChar(coeffs cf, void * infoStruct)
1674{ return TRUE; }
1675#endif
1676#else
1677BOOLEAN flintQrat_InitChar(coeffs cf, void * infoStruct)
1678{ return TRUE; }
1679#endif
Note: See TracBrowser for help on using the repository browser.