My Project
Loading...
Searching...
No Matches
polys.cc
Go to the documentation of this file.
1#include "kernel/mod2.h"
2
3#include "misc/options.h"
4
5#include "polys.h"
6#include "kernel/ideals.h"
7#include "kernel/ideals.h"
8#include "polys/clapsing.h"
9#include "polys/clapconv.h"
10
11/// Widely used global variable which specifies the current polynomial ring for Singular interpreter and legacy implementations.
12/// @Note: one should avoid using it in newer designs, for example due to possible problems in parallelization with threads.
15
16void rChangeCurrRing(ring r)
17{
18 if (currRing!=NULL)
20 //------------ set global ring vars --------------------------------
21 currRing = r;
22 if( r != NULL )
23 {
24 rTest(r);
25 //------------ global variables related to coefficients ------------
26 assume( r->cf!= NULL );
27 nSetChar(r->cf);
28 //------------ global variables related to polys
29 p_SetGlobals(r); // also setting TEST_RINGDEP_OPTS
30 //------------ global variables related to factory -----------------
31 }
32}
33
34poly p_Divide(poly p, poly q, const ring r)
35{
36 assume(q!=NULL);
37 if (q==NULL)
38 {
39 WerrorS("div. by 0");
40 return NULL;
41 }
42 if (p==NULL)
43 {
44 p_Delete(&q,r);
45 return NULL;
46 }
47 if ((pNext(q)!=NULL)||rIsPluralRing(r))
48 { /* This means that q != 0 consists of at least two terms*/
49 if(p_GetComp(p,r)==0)
50 {
51 if (!rIsNCRing(r))
52 {
53 if ((rFieldType(r)==n_Q)
54 ||(rFieldType(r)==n_Zp)
55 ||(rFieldType(r)==n_algExt))
56 {
57 poly res=singclap_pdivide(p, q, r);
58 p_Delete(&p,r);
59 p_Delete(&q,r);
60 return res;
61 }
62 if((rFieldType(r)==n_transExt)
63 &&(convSingTrP(p,r))
64 &&(convSingTrP(q,r)))
65 {
66 poly res=singclap_pdivide(p, q, r);
67 p_Delete(&p,r);
68 p_Delete(&q,r);
69 return res;
70 }
71 }
72 // generic division for poly
73 {
74 ideal vi=idInit(1,1); vi->m[0]=q;
75 ideal ui=idInit(1,1); ui->m[0]=p;
76 ideal R; matrix U;
77 ring save_ring=currRing;
78 if (r!=currRing) rChangeCurrRing(r);
79 BITSET save_opt;
80 SI_SAVE_OPT1(save_opt);
82 ideal m = idLift(vi,ui,&R, FALSE,TRUE,TRUE,&U);
83 SI_RESTORE_OPT1(save_opt);
84 if (r!=save_ring) rChangeCurrRing(save_ring);
86 p=MATELEM(T,1,1); MATELEM(T,1,1)=NULL;
87 id_Delete((ideal *)&T,r);
88 id_Delete((ideal *)&U,r);
89 id_Delete(&R,r);
90 //vi->m[0]=NULL; ui->m[0]=NULL;
91 id_Delete(&vi,r);
92 id_Delete(&ui,r);
93 return p;
94 }
95 }
96 else
97 {
98 int comps=p_MaxComp(p,r);
99 ideal I=idInit(comps,1);
100 poly h;
101 int i;
102 // conversion to a list of polys:
103 while (p!=NULL)
104 {
105 i=p_GetComp(p,r)-1;
106 h=pNext(p);
107 pNext(p)=NULL;
108 p_SetComp(p,0,r);
109 I->m[i]=p_Add_q(I->m[i],p,r);
110 p=h;
111 }
112 // division and conversion to vector:
113 h=NULL;
114 p=NULL;
115 for(i=comps-1;i>=0;i--)
116 {
117 if (I->m[i]!=NULL)
118 {
119 if((rFieldType(r)==n_transExt)
120 &&(convSingTrP(I->m[i],r))
121 &&(convSingTrP(q,r))
122 &&(!rIsNCRing(r)))
123 {
124 h=singclap_pdivide(I->m[i],q,r);
125 }
126 else if (((rFieldType(r)==n_Q)
127 ||(rFieldType(r)==n_Zp)
128 ||(rFieldType(r)==n_algExt))
129 &&(!rIsNCRing(r)))
130 h=singclap_pdivide(I->m[i],q,r);
131 else
132 {
133 ideal vi=idInit(1,1); vi->m[0]=q;
134 ideal ui=idInit(1,1); ui->m[0]=I->m[i];
135 ideal R; matrix U;
136 ring save_ring=currRing;
137 if (r!=currRing) rChangeCurrRing(r);
138 BITSET save_opt;
139 SI_SAVE_OPT1(save_opt);
140 si_opt_1 &= ~(Sy_bit(OPT_PROT));
141 ideal m = idLift(vi,ui,&R, FALSE,TRUE,TRUE,&U);
142 SI_RESTORE_OPT1(save_opt);
143 if (r!=save_ring) rChangeCurrRing(save_ring);
144 if (idIs0(R))
145 {
147 p=MATELEM(T,1,1); MATELEM(T,1,1)=NULL;
148 id_Delete((ideal *)&T,r);
149 }
150 else p=NULL;
151 id_Delete((ideal*)&U,r);
152 id_Delete(&R,r);
153 vi->m[0]=NULL; ui->m[0]=NULL;
154 id_Delete(&vi,r);
155 id_Delete(&ui,r);
156 }
157 p_SetCompP(h,i+1,r);
158 p=p_Add_q(p,h,r);
159 }
160 }
161 id_Delete(&I,r);
162 p_Delete(&q,r);
163 return p;
164 }
165 }
166 else
167 { /* This means that q != 0 consists of just one term, or LetterPlace */
168 if (pNext(q)!=NULL)
169 {
170 WerrorS("division over a coefficient domain only implemented for terms");
171 return NULL;
172 }
173 return p_DivideM(p,q,r);
174 }
175 return NULL;
176}
177
178poly pp_Divide(poly p, poly q, const ring r)
179{
180 if (q==NULL)
181 {
182 WerrorS("div. by 0");
183 return NULL;
184 }
185 if (p==NULL)
186 {
187 return NULL;
188 }
189 if ((pNext(q)!=NULL)||rIsPluralRing(r))
190 { /* This means that q != 0 consists of at least two terms*/
191 if(p_GetComp(p,r)==0)
192 {
193 if (!rIsNCRing(r))
194 {
195 if ((rFieldType(r)==n_Q)
196 ||(rFieldType(r)==n_Zp)
197 ||(rFieldType(r)==n_algExt))
198 {
199 poly res=singclap_pdivide(p, q, r);
200 return res;
201 }
202 if((rFieldType(r)==n_transExt)
203 &&(convSingTrP(p,r))
204 &&(convSingTrP(q,r)))
205 {
206 poly res=singclap_pdivide(p, q, r);
207 return res;
208 }
209 }
210 // generic division for poly
211 {
212 ideal vi=idInit(1,1); vi->m[0]=p_Copy(q,r);
213 ideal ui=idInit(1,1); ui->m[0]=p_Copy(p,r);
214 ideal R; matrix U;
215 ring save_ring=currRing;
216 if (r!=currRing) rChangeCurrRing(r);
217 BITSET save_opt;
218 SI_SAVE_OPT1(save_opt);
219 si_opt_1 &= ~(Sy_bit(OPT_PROT));
220 ideal m = idLift(vi,ui,&R, FALSE,TRUE,TRUE,&U);
221 SI_RESTORE_OPT1(save_opt);
222 if (r!=save_ring) rChangeCurrRing(save_ring);
224 p=MATELEM(T,1,1); MATELEM(T,1,1)=NULL;
225 id_Delete((ideal *)&T,r);
226 id_Delete((ideal *)&U,r);
227 id_Delete(&R,r);
228 //vi->m[0]=NULL; ui->m[0]=NULL;
229 id_Delete(&vi,r);
230 id_Delete(&ui,r);
231 return p;
232 }
233 }
234 else
235 {
236 p=p_Copy(p,r);
237 int comps=p_MaxComp(p,r);
238 ideal I=idInit(comps,1);
239 poly h;
240 int i;
241 // conversion to a list of polys:
242 while (p!=NULL)
243 {
244 i=p_GetComp(p,r)-1;
245 h=pNext(p);
246 pNext(p)=NULL;
247 p_SetComp(p,0,r);
248 I->m[i]=p_Add_q(I->m[i],p,r);
249 p=h;
250 }
251 // division and conversion to vector:
252 h=NULL;
253 p=NULL;
254 q=p_Copy(q,r);
255 for(i=comps-1;i>=0;i--)
256 {
257 if (I->m[i]!=NULL)
258 {
259 if(!rIsNCRing(r))
260 {
261 if((rFieldType(r)==n_transExt)
262 &&(convSingTrP(I->m[i],r))
263 &&(convSingTrP(q,r)))
264 {
265 h=singclap_pdivide(I->m[i],q,r);
266 }
267 else if ((rFieldType(r)==n_Q)
268 ||(rFieldType(r)==n_Zp))
269 h=singclap_pdivide(I->m[i],q,r);
270 else
271 {
272 ideal vi=idInit(1,1); vi->m[0]=q;
273 ideal ui=idInit(1,1); ui->m[0]=I->m[i];
274 ideal R; matrix U;
275 ring save_ring=currRing;
276 if (r!=currRing) rChangeCurrRing(r);
277 BITSET save_opt;
278 SI_SAVE_OPT1(save_opt);
279 si_opt_1 &= ~(Sy_bit(OPT_PROT));
280 ideal m = idLift(vi,ui,&R, FALSE,TRUE,TRUE,&U);
281 SI_RESTORE_OPT1(save_opt);
282 if (r!=save_ring) rChangeCurrRing(save_ring);
283 if (idIs0(R))
284 {
286 p=MATELEM(T,1,1); MATELEM(T,1,1)=NULL;
287 id_Delete((ideal *)&T,r);
288 }
289 id_Delete((ideal*)&U,r);
290 id_Delete(&R,r);
291 vi->m[0]=NULL; ui->m[0]=NULL;
292 id_Delete(&vi,r);
293 id_Delete(&ui,r);
294 }
295 }
296 else
297 {
298 ideal vi=idInit(1,1); vi->m[0]=q;
299 ideal ui=idInit(1,1); ui->m[0]=I->m[i];
300 ideal R; matrix U;
301 ring save_ring=currRing;
302 if (r!=currRing) rChangeCurrRing(r);
303 BITSET save_opt;
304 SI_SAVE_OPT1(save_opt);
305 si_opt_1 &= ~(Sy_bit(OPT_PROT));
306 ideal m = idLift(vi,ui,&R, FALSE,TRUE,TRUE,&U);
307 SI_RESTORE_OPT1(save_opt);
308 if (r!=save_ring) rChangeCurrRing(save_ring);
309 if (idIs0(R))
310 {
312 p=MATELEM(T,1,1); MATELEM(T,1,1)=NULL;
313 id_Delete((ideal *)&T,r);
314 }
315 id_Delete((ideal*)&U,r);
316 id_Delete(&R,r);
317 vi->m[0]=NULL; ui->m[0]=NULL;
318 id_Delete(&vi,r);
319 id_Delete(&ui,r);
320 }
321 p_SetCompP(h,i+1,r);
322 p=p_Add_q(p,h,r);
323 }
324 }
325 id_Delete(&I,r);
326 p_Delete(&q,r);
327 return p;
328 }
329 }
330 else
331 { /* This means that q != 0 consists of just one term,
332 or that r is over a coefficient ring. */
333 if (pNext(q)!=NULL)
334 {
335 WerrorS("division over a coefficient domain only implemented for terms");
336 return NULL;
337 }
338 return pp_DivideM(p,q,r);
339 }
340 return NULL;
341}
342
343poly p_DivRem(poly p, poly q, poly &rest, const ring r)
344{
345 assume(q!=NULL);
346 rest=NULL;
347 if (q==NULL)
348 {
349 WerrorS("div. by 0");
350 return NULL;
351 }
352 if (p==NULL)
353 {
354 p_Delete(&q,r);
355 return NULL;
356 }
357 if(p_GetComp(p,r)==0)
358 {
359 if((rFieldType(r)==n_transExt)
360 &&(convSingTrP(p,r))
361 &&(convSingTrP(q,r))
362 &&(!rIsNCRing(r)))
363 {
364 poly res=singclap_pdivide(p, q, r);
365 rest=singclap_pmod(p,q,r);
366 p_Delete(&p,r);
367 p_Delete(&q,r);
368 return res;
369 }
370 else if ((r->cf->convSingNFactoryN!=ndConvSingNFactoryN)
371 &&(!rField_is_Ring(r))
372 &&(!rIsNCRing(r)))
373 {
374 poly res=singclap_pdivide(p, q, r);
375 rest=singclap_pmod(p,q,r);
376 p_Delete(&p,r);
377 p_Delete(&q,r);
378 return res;
379 }
380 else
381 {
382 ideal vi=idInit(1,1); vi->m[0]=q;
383 ideal ui=idInit(1,1); ui->m[0]=p;
384 ideal R; matrix U;
385 ring save_ring=currRing;
386 if (r!=currRing) rChangeCurrRing(r);
387 BITSET save_opt;
388 SI_SAVE_OPT1(save_opt);
389 si_opt_1 &= ~(Sy_bit(OPT_PROT));
390 ideal m = idLift(vi,ui,&R, FALSE,TRUE,TRUE,&U);
391 SI_RESTORE_OPT1(save_opt);
392 if (r!=save_ring) rChangeCurrRing(save_ring);
393 p=m->m[0]; m->m[0]=NULL;
394 id_Delete(&m,r);
395 p_SetCompP(p,0,r);
396 rest=R->m[0]; R->m[0]=NULL;
397 id_Delete(&R,r);
398 p_SetCompP(rest,0,r);
399 id_Delete((ideal *)&U,r);
400 //vi->m[0]=NULL; ui->m[0]=NULL;
401 id_Delete(&vi,r);
402 id_Delete(&ui,r);
403 return p;
404 }
405 }
406 return NULL;
407}
408
409poly singclap_gcd ( poly f, poly g, const ring r )
410{
411 poly res=NULL;
412
413 if (f!=NULL)
414 {
415 //if (r->cf->has_simple_Inverse) p_Norm(f,r);
416 if (rField_is_Zp(r)) p_Norm(f,r);
417 else if (!rField_is_Ring(r)) p_Cleardenom(f, r);
418 }
419 if (g!=NULL)
420 {
421 //if (r->cf->has_simple_Inverse) p_Norm(g,r);
422 if (rField_is_Zp(r)) p_Norm(g,r);
423 else if (!rField_is_Ring(r)) p_Cleardenom(g, r);
424 }
425 else return f; // g==0 => gcd=f (but do a p_Cleardenom/pNorm)
426 if (f==NULL) return g; // f==0 => gcd=g (but do a p_Cleardenom/pNorm)
427 if(!rField_is_Ring(r)
428 && (p_IsConstant(f,r)
429 ||p_IsConstant(g,r)))
430 {
431 res=p_One(r);
432 }
433 else if (r->cf->convSingNFactoryN!=ndConvSingNFactoryN)
434 {
436 }
437 else
438 {
439 ideal I=idInit(2,1);
440 I->m[0]=f;
441 I->m[1]=p_Copy(g,r);
442 intvec *w=NULL;
443 ring save_ring=currRing;
444 if (r!=currRing) rChangeCurrRing(r);
445 BITSET save_opt;
446 SI_SAVE_OPT1(save_opt);
447 si_opt_1 &= ~(Sy_bit(OPT_PROT));
448 ideal S1=idSyzygies(I,testHomog,&w);
449 if (w!=NULL) delete w;
450 // expect S1->m[0]=(-g/gcd,f/gcd)
451 if (IDELEMS(S1)!=1) WarnS("error in syzygy computation for GCD");
452 int lp;
453 p_TakeOutComp(&S1->m[0],1,&res,&lp,r);
454 p_Delete(&S1->m[0],r);
455 // GCD is g divided iby (-g/gcd):
456 res=p_Divide(g,res,r);
457 // restore, r, opt:
458 SI_RESTORE_OPT1(save_opt);
459 if (r!=save_ring) rChangeCurrRing(save_ring);
460 // clean the result
462 if (nCoeff_is_Ring(r->cf)) p_Content(res,r);
463 return res;
464 }
465 p_Delete(&f, r);
466 p_Delete(&g, r);
467 return res;
468}
#define BITSET
Definition auxiliary.h:85
#define TRUE
Definition auxiliary.h:101
#define FALSE
Definition auxiliary.h:97
int m
Definition cfEzgcd.cc:128
int i
Definition cfEzgcd.cc:132
int p
Definition cfModGcd.cc:4086
g
Definition cfModGcd.cc:4098
FILE * f
Definition checklibs.c:9
BOOLEAN convSingTrP(poly p, const ring r)
Definition clapconv.cc:375
poly singclap_pmod(poly f, poly g, const ring r)
Definition clapsing.cc:732
poly singclap_pdivide(poly f, poly g, const ring r)
Definition clapsing.cc:649
poly singclap_gcd_r(poly f, poly g, const ring r)
Definition clapsing.cc:70
@ n_Q
rational (GMP) numbers
Definition coeffs.h:30
@ n_algExt
used for all algebraic extensions, i.e., the top-most extension in an extension tower is algebraic
Definition coeffs.h:35
@ n_Zp
\F{p < 2^31}
Definition coeffs.h:29
@ n_transExt
used for all transcendental extensions, i.e., the top-most extension in an extension tower is transce...
Definition coeffs.h:38
static FORCE_INLINE void nSetChar(const coeffs r)
initialisations after each ring change
Definition coeffs.h:446
static FORCE_INLINE BOOLEAN nCoeff_is_Ring(const coeffs r)
Definition coeffs.h:732
#define WarnS
Definition emacs.cc:78
CanonicalForm res
Definition facAbsFact.cc:60
const CanonicalForm & w
Definition facAbsFact.cc:51
void WerrorS(const char *s)
Definition feFopen.cc:24
#define VAR
Definition globaldefs.h:5
ideal idSyzygies(ideal h1, tHomog h, intvec **w, BOOLEAN setSyzComp, BOOLEAN setRegularity, int *deg, GbVariant alg)
Definition ideals.cc:834
ideal idLift(ideal mod, ideal submod, ideal *rest, BOOLEAN goodShape, BOOLEAN isSB, BOOLEAN divide, matrix *unit, GbVariant alg)
represents the generators of submod in terms of the generators of mod (Matrix(SM)*U-Matrix(rest)) = M...
Definition ideals.cc:1109
BOOLEAN idIs0(ideal h)
returns true if h is the zero ideal
STATIC_VAR jList * T
Definition janet.cc:30
STATIC_VAR Poly * h
Definition janet.cc:971
void p_TakeOutComp(poly *p, long comp, poly *q, int *lq, const ring r)
Definition p_polys.cc:3620
#define MATELEM(mat, i, j)
1-based access to matrix
Definition matpol.h:29
ip_smatrix * matrix
Definition matpol.h:43
#define assume(x)
Definition mod2.h:389
#define p_GetComp(p, r)
Definition monomials.h:64
#define pNext(p)
Definition monomials.h:36
The main handler for Singular numbers which are suitable for Singular polynomials.
CanonicalForm ndConvSingNFactoryN(number, BOOLEAN, const coeffs)
Definition numbers.cc:313
#define NULL
Definition omList.c:12
VAR unsigned si_opt_1
Definition options.c:5
#define OPT_PROT
Definition options.h:76
#define SI_SAVE_OPT1(A)
Definition options.h:21
#define SI_RESTORE_OPT1(A)
Definition options.h:24
#define Sy_bit(x)
Definition options.h:31
#define TEST_RINGDEP_OPTS
Definition options.h:101
poly p_DivideM(poly a, poly b, const ring r)
Definition p_polys.cc:1582
void p_Content(poly ph, const ring r)
Definition p_polys.cc:2343
poly pp_DivideM(poly a, poly b, const ring r)
Definition p_polys.cc:1637
void p_Norm(poly p1, const ring r)
Definition p_polys.cc:3844
poly p_Cleardenom(poly p, const ring r)
Definition p_polys.cc:2893
poly p_One(const ring r)
Definition p_polys.cc:1314
static poly p_Add_q(poly p, poly q, const ring r)
Definition p_polys.h:938
static void p_SetCompP(poly p, int i, ring r)
Definition p_polys.h:256
static unsigned long p_SetComp(poly p, unsigned long c, ring r)
Definition p_polys.h:249
static BOOLEAN p_IsConstant(const poly p, const ring r)
Definition p_polys.h:1985
static long p_MaxComp(poly p, ring lmRing, ring tailRing)
Definition p_polys.h:294
static void p_Delete(poly *p, const ring r)
Definition p_polys.h:903
static poly p_Copy(poly p, const ring r)
returns a copy of p
Definition p_polys.h:848
void rChangeCurrRing(ring r)
Definition polys.cc:16
VAR coeffs coeffs_BIGINT
Definition polys.cc:14
poly p_Divide(poly p, poly q, const ring r)
polynomial division a/b, ignoring the rest via singclap_pdivide resp. idLift destroys a,...
Definition polys.cc:34
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition polys.cc:13
poly pp_Divide(poly p, poly q, const ring r)
polynomial division a/b, ignoring the rest via singclap_pdivide resp. idLift does not destroy a,...
Definition polys.cc:178
poly singclap_gcd(poly f, poly g, const ring r)
polynomial gcd via singclap_gcd_r resp. idSyzygies destroys f and g
Definition polys.cc:409
poly p_DivRem(poly p, poly q, poly &rest, const ring r)
Definition polys.cc:343
Compatibility layer for legacy polynomial operations (over currRing).
void p_SetGlobals(const ring r, BOOLEAN complete)
set all properties of a new ring - also called by rComplete
Definition ring.cc:3494
static BOOLEAN rField_is_Zp(const ring r)
Definition ring.h:506
static n_coeffType rFieldType(const ring r)
the type of the coefficient filed of r (n_Zp, n_Q, etc)
Definition ring.h:567
static BOOLEAN rIsPluralRing(const ring r)
we must always have this test!
Definition ring.h:406
static BOOLEAN rIsNCRing(const ring r)
Definition ring.h:427
#define rTest(r)
Definition ring.h:799
#define rField_is_Ring(R)
Definition ring.h:491
ideal idInit(int idsize, int rank)
initialise an ideal / module
void id_Delete(ideal *h, ring r)
deletes an ideal/module/matrix
matrix id_Module2formatedMatrix(ideal mod, int rows, int cols, const ring R)
#define IDELEMS(i)
#define R
Definition sirandom.c:27
@ testHomog
Definition structs.h:34