source: git/libpolys/polys/flint_mpoly.cc @ 4abcd2

fieker-DuValspielwiese
Last change on this file since 4abcd2 was 4abcd2, checked in by Hans Schoenemann <hannes@…>, 5 years ago
add: gcd via flint (over QQ,ZZ/p)
  • Property mode set to 100644
File size: 5.9 KB
Line 
1// emacs edit mode for this file is -*- C++ -*-
2/****************************************
3*  Computer Algebra System SINGULAR     *
4****************************************/
5/*
6* ABSTRACT: flint mpoly
7*/
8
9#include "misc/auxiliary.h"
10#include "flintconv.h"
11#include "flint_mpoly.h"
12
13#ifdef HAVE_FLINT
14#if __FLINT_RELEASE >= 20503
15#include "coeffs/coeffs.h"
16#include "coeffs/longrat.h"
17#include "polys/monomials/p_polys.h"
18
19BOOLEAN convSingRFlintR(fmpq_mpoly_ctx_t ctx, const ring r)
20{
21  if (rRing_ord_pure_dp(r))
22  {
23    fmpq_mpoly_ctx_init(ctx,r->N,ORD_DEGREVLEX);
24    return FALSE;
25  }
26  else if (rRing_ord_pure_Dp(r))
27  {
28    fmpq_mpoly_ctx_init(ctx,r->N,ORD_DEGLEX);
29    return FALSE;
30  }
31  else if (rRing_ord_pure_lp(r))
32  {
33    fmpq_mpoly_ctx_init(ctx,r->N,ORD_LEX);
34    return FALSE;
35  }
36  return TRUE;
37}
38
39void convSingPFlintMP(fmpq_mpoly_t res, fmpq_mpoly_ctx_t ctx, poly p, int lp, const ring r)
40{
41  int bits=SI_LOG2(r->bitmask);
42  fmpq_mpoly_init3(res,lp,bits,ctx);
43  fmpq_mpoly_resize(res,lp,ctx);
44  ulong* exp=(ulong*)omAlloc((r->N+1)*sizeof(ulong));
45  int i=0;
46  while(p!=NULL)
47  {
48    number n=pGetCoeff(p);
49    fmpq_t c;
50    convSingNFlintN_QQ(c,n);
51    #if SIZEOF_LONG==8
52    p_GetExpVL(p,(int64*)exp,r);
53    fmpq_mpoly_set_term_exp_ui(res,i,exp,ctx);
54    #else
55    p_GetExpV(p,(int*)exp,r);
56    fmpq_mpoly_set_term_exp_ui(res,i,&(exp[1]),ctx);
57    #endif
58    fmpq_mpoly_set_term_coeff_fmpq(res,i,c,ctx);
59    fmpq_clear(c);
60    i++;
61    pIter(p);
62  }
63  //fmpq_mpoly_print_pretty(res,r->names,ctx);PrintLn();
64  omFreeSize(exp,(r->N+1)*sizeof(ulong));
65}
66
67poly convFlintMPSingP(fmpq_mpoly_t f, fmpq_mpoly_ctx_t ctx, const ring r)
68{
69  int d=fmpq_mpoly_length(f,ctx)-1;
70  poly p=NULL;
71  ulong* exp=(ulong*)omAlloc0((r->N+1)*sizeof(ulong));
72  for(int i=d; i>=0; i--)
73  {
74    fmpq_t c;
75    fmpq_mpoly_get_term_coeff_fmpq(c,f,i,ctx);
76    poly pp=p_Init(r);
77    #if SIZEOF_LONG==8
78    fmpq_mpoly_get_term_exp_ui(exp,f,i,ctx);
79    p_SetExpVL(pp,(int64*)exp,r);
80    #else
81    fmpq_mpoly_get_term_exp_ui(&(exp[1]),f,i,ctx);
82    p_SetExpV(pp,(int*)exp,r);
83    #endif
84    p_Setm(pp,r);
85    number n=convFlintNSingN_QQ(c,r->cf);
86    //fmpq_clear(c); // LEAK?
87    pSetCoeff0(pp,n);
88    pNext(pp)=p;
89    p=pp;
90  }
91  p_Test(p,r);
92  return p;
93}
94
95poly Flint_Mult_MP(poly p,int lp, poly q, int lq, fmpq_mpoly_ctx_t ctx, const ring r)
96{
97  fmpq_mpoly_t pp,qq,res;
98  convSingPFlintMP(pp,ctx,p,lp,r);
99  convSingPFlintMP(qq,ctx,q,lq,r);
100  int bits=SI_LOG2(r->bitmask);
101  fmpq_mpoly_init3(res,lp*lq,bits,ctx);
102  fmpq_mpoly_mul(res,pp,qq,ctx);
103  poly pres=convFlintMPSingP(res,ctx,r);
104  fmpq_mpoly_clear(res,ctx);
105  fmpq_mpoly_clear(pp,ctx);
106  fmpq_mpoly_clear(qq,ctx);
107  fmpq_mpoly_ctx_clear(ctx);
108  p_Test(pres,r);
109  return pres;
110}
111BOOLEAN convSingRFlintR(nmod_mpoly_ctx_t ctx, const ring r)
112{
113  if (rRing_ord_pure_dp(r))
114  {
115    nmod_mpoly_ctx_init(ctx,r->N,ORD_DEGREVLEX,r->cf->ch);
116    return FALSE;
117  }
118  else if (rRing_ord_pure_Dp(r))
119  {
120    nmod_mpoly_ctx_init(ctx,r->N,ORD_DEGLEX,r->cf->ch);
121    return FALSE;
122  }
123  else if (rRing_ord_pure_lp(r))
124  {
125    nmod_mpoly_ctx_init(ctx,r->N,ORD_LEX,r->cf->ch);
126    return FALSE;
127  }
128  return TRUE;
129}
130
131poly convFlintMPSingP(nmod_mpoly_t f, nmod_mpoly_ctx_t ctx, const ring r)
132{
133  int d=nmod_mpoly_length(f,ctx)-1;
134  poly p=NULL;
135  ulong* exp=(ulong*)omAlloc0((r->N+1)*sizeof(ulong));
136  for(int i=d; i>=0; i--)
137  {
138    ulong c=nmod_mpoly_get_term_coeff_ui(f,i,ctx);
139    poly pp=p_Init(r);
140    #if SIZEOF_LONG==8
141    nmod_mpoly_get_term_exp_ui(exp,f,i,ctx);
142    p_SetExpVL(pp,(int64*)exp,r);
143    #else
144    nmod_mpoly_get_term_exp_ui(&(exp[1]),f,i,ctx);
145    p_SetExpV(pp,(int*)exp,r);
146    #endif
147    p_Setm(pp,r);
148    pSetCoeff0(pp,(number)c);
149    pNext(pp)=p;
150    p=pp;
151  }
152  p_Test(p,r);
153  return p;
154}
155
156void convSingPFlintMP(nmod_mpoly_t res, nmod_mpoly_ctx_t ctx, poly p, int lp,const ring r)
157{
158  int bits=SI_LOG2(r->bitmask);
159  nmod_mpoly_init3(res,lp,bits,ctx);
160  nmod_mpoly_resize(res,lp,ctx);
161  ulong* exp=(ulong*)omAlloc((r->N+1)*sizeof(ulong));
162  int i=0;
163  while(p!=NULL)
164  {
165    number n=pGetCoeff(p);
166    #if SIZEOF_LONG==8
167    p_GetExpVL(p,(int64*)exp,r);
168    nmod_mpoly_set_term_exp_ui(res,i,exp,ctx);
169    #else
170    p_GetExpV(p,(int*)exp,r);
171    nmod_mpoly_set_term_exp_ui(res,i,&(exp[1]),ctx);
172    #endif
173    nmod_mpoly_set_term_coeff_ui(res,i,(ulong)n,ctx);
174    i++;
175    pIter(p);
176  }
177  //nmod_mpoly_print_pretty(res,r->names,ctx);PrintLn();
178  omFreeSize(exp,(r->N+1)*sizeof(ulong));
179}
180
181poly Flint_Mult_MP(poly p,int lp, poly q, int lq, nmod_mpoly_ctx_t ctx, const ring r)
182{
183  nmod_mpoly_t pp,qq,res;
184  convSingPFlintMP(pp,ctx,p,lp,r);
185  convSingPFlintMP(qq,ctx,q,lq,r);
186  int bits=SI_LOG2(r->bitmask);
187  nmod_mpoly_init3(res,lp*lq,bits,ctx);
188  nmod_mpoly_mul(res,pp,qq,ctx);
189  poly pres=convFlintMPSingP(res,ctx,r);
190  nmod_mpoly_clear(res,ctx);
191  nmod_mpoly_clear(pp,ctx);
192  nmod_mpoly_clear(qq,ctx);
193  nmod_mpoly_ctx_clear(ctx);
194  p_Test(pres,r);
195  return pres;
196}
197poly Flint_GCD_MP(poly p,int lp,poly q,int lq,nmod_mpoly_ctx_t ctx,const ring r)
198{
199  nmod_mpoly_t pp,qq,res;
200  convSingPFlintMP(pp,ctx,p,lp,r);
201  convSingPFlintMP(qq,ctx,q,lq,r);
202  int bits=SI_LOG2(r->bitmask);
203  nmod_mpoly_init3(res,lp*lq,bits,ctx);
204  int ok=nmod_mpoly_gcd(res,pp,qq,ctx);
205  poly pres;
206  if (ok)
207  {
208    pres=convFlintMPSingP(res,ctx,r);
209    p_Test(pres,r);
210  }
211  else
212  {
213    pres=p_One(r);
214  }
215  nmod_mpoly_clear(res,ctx);
216  nmod_mpoly_clear(pp,ctx);
217  nmod_mpoly_clear(qq,ctx);
218  nmod_mpoly_ctx_clear(ctx);
219  return pres;
220}
221
222poly Flint_GCD_MP(poly p,int lp,poly q,int lq,fmpq_mpoly_ctx_t ctx,const ring r)
223{
224  fmpq_mpoly_t pp,qq,res;
225  convSingPFlintMP(pp,ctx,p,lp,r);
226  convSingPFlintMP(qq,ctx,q,lq,r);
227  int bits=SI_LOG2(r->bitmask);
228  fmpq_mpoly_init3(res,lp*lq,bits,ctx);
229  int ok=fmpq_mpoly_gcd(res,pp,qq,ctx);
230  poly pres;
231  if (ok)
232  {
233    pres=convFlintMPSingP(res,ctx,r);
234    p_Test(pres,r);
235  }
236  else
237  {
238    pres=p_One(r);
239  }
240  fmpq_mpoly_clear(res,ctx);
241  fmpq_mpoly_clear(pp,ctx);
242  fmpq_mpoly_clear(qq,ctx);
243  fmpq_mpoly_ctx_clear(ctx);
244  return pres;
245}
246#endif
247#endif
Note: See TracBrowser for help on using the repository browser.