My Project
clapsing.cc
Go to the documentation of this file.
1// emacs edit mode for this file is -*- C++ -*-
2/****************************************
3* Computer Algebra System SINGULAR *
4****************************************/
5/*
6* ABSTRACT: interface between Singular and factory
7*/
8
9//#define FACTORIZE2_DEBUG
10
11#include "misc/auxiliary.h"
12#include "clapsing.h"
13
14#include "factory/factory.h"
15#include "factory/cf_roots.h"
16
17#include "coeffs/numbers.h"
18#include "coeffs/coeffs.h"
19#include "coeffs/bigintmat.h"
20
21#include "monomials/ring.h"
22#include "simpleideals.h"
23#include "polys/flintconv.h"
24#include "polys/flint_mpoly.h"
25
26#ifdef HAVE_NTL
27#include <NTL/config.h>
28#ifdef NTL_STD_CXX
29#ifdef NOSTREAMIO
30# ifdef HAVE_IOSTREAM
31# include <iostream>
32# define OSTREAM std::ostream
33# define ISTREAM std::istream
34# elif defined(HAVE_IOSTREAM_H)
35# include <iostream.h>
36# define OSTREAM ostream
37# define ISTREAM istream
38# endif
39#endif /* ! NOSTREAMIO */
40#endif
41#include <NTL/mat_ZZ.h>
42#include <NTL/mat_lzz_p.h>
43
44#ifdef NTL_CLIENT // in <NTL/tools.h>: using of name space NTL
45NTL_CLIENT
46#endif
47
48#endif
49
50//#include "polys.h"
51#define TRANSEXT_PRIVATES
52
53#include "ext_fields/transext.h"
54
55
56#include "clapconv.h"
57
58#include "polys/monomials/p_polys.h"
59#include "polys/monomials/ring.h"
60#include "polys/simpleideals.h"
61#include "misc/intvec.h"
62#include "polys/matpol.h"
63#include "coeffs/bigintmat.h"
64
65
66void out_cf(const char *s1,const CanonicalForm &f,const char *s2);
67
68poly singclap_gcd_r ( poly f, poly g, const ring r )
69{
70 poly res=NULL;
71
72 assume(f!=NULL);
73 assume(g!=NULL);
74
75 if(pNext(f)==NULL)
76 {
77 return p_GcdMon(f,g,r);
78 }
79 else if(pNext(g)==NULL)
80 {
81 return p_GcdMon(g,f,r);
82 }
83 #ifdef HAVE_FLINT
84 #if __FLINT_RELEASE >= 20503
85 if (rField_is_Zp(r) && (r->cf->ch>10))
86 {
87 nmod_mpoly_ctx_t ctx;
88 if (!convSingRFlintR(ctx,r))
89 {
90 // leading coef. 1
91 return Flint_GCD_MP(f,pLength(f),g,pLength(g),ctx,r);
92 }
93 }
94 else
95 if (rField_is_Q(r))
96 {
97 fmpq_mpoly_ctx_t ctx;
98 if (!convSingRFlintR(ctx,r))
99 {
100 // leading coef. positive, all coeffs in Z
101 poly res=Flint_GCD_MP(f,pLength(f),g,pLength(g),ctx,r);
102 res=p_Cleardenom(res,r);
103 return res;
104 }
105 }
106 else
107 if (rField_is_Z(r))
108 {
109 fmpz_mpoly_ctx_t ctx;
110 if (!convSingRFlintR(ctx,r))
111 {
112 // leading coef. positive, all coeffs in Z
113 poly res=Flint_GCD_MP(f,pLength(f),g,pLength(g),ctx,r);
114 return res;
115 }
116 }
117 #endif
118 #endif
119 Off(SW_RATIONAL);
120 if (rField_is_Q(r) || rField_is_Zp(r) || rField_is_Z(r)
121 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
122 {
123 setCharacteristic( rInternalChar(r) );
124 CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g, r ) );
125 res=convFactoryPSingP( gcd( F, G ) , r);
126 if ( rField_is_Zp(r))
127 p_Norm(res,r); // leading coef. 1
128 else if (rField_is_Q(r) && (!n_GreaterZero(pGetCoeff(res),r->cf)))
129 res = p_Neg(res,r); // leading coef. positive, all coeffs in Z
130 }
131 // and over Q(a) / Fp(a)
132 else if ( r->cf->extRing!=NULL )
133 {
134 if ( rField_is_Q_a(r)) setCharacteristic( 0 );
135 else setCharacteristic( rInternalChar(r) );
136 if (r->cf->extRing->qideal!=NULL)
137 {
138 bool b1=isOn(SW_USE_QGCD);
139 if ( rField_is_Q_a(r) ) On(SW_USE_QGCD);
140 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
141 r->cf->extRing);
142 Variable a=rootOf(mipo);
143 CanonicalForm F( convSingAPFactoryAP( f,a,r ) ),
144 G( convSingAPFactoryAP( g,a,r ) );
145 res= convFactoryAPSingAP( gcd( F, G ),r );
146 prune (a);
147 if (!b1) Off(SW_USE_QGCD);
148 if ( rField_is_Zp_a(r)) p_Norm(res,r); // leading coef. 1
149 }
150 else
151 {
152 convSingTrP(f,r);
153 convSingTrP(g,r);
154 CanonicalForm F( convSingTrPFactoryP( f,r ) ), G( convSingTrPFactoryP( g,r ) );
155 res= convFactoryPSingTrP( gcd( F, G ),r );
156 }
157 }
158 else if (r->cf->convSingNFactoryN==ndConvSingNFactoryN)
159 WerrorS( feNotImplemented );
160 else
161 { // handle user type coeffs:
162 setCharacteristic( rInternalChar(r) );
163 CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g, r ) );
164 res=convFactoryPSingP( gcd( F, G ) , r);
165 }
166 Off(SW_RATIONAL);
167 return res;
168}
169
170poly singclap_gcd_and_divide ( poly& f, poly& g, const ring r)
171{
172 poly res=NULL;
173
174 if (g == NULL)
175 {
176 res= f;
177 f=p_One (r);
178 return res;
179 }
180 if (f==NULL)
181 {
182 res= g;
183 g=p_One (r);
184 return res;
185 }
186 if (pNext(g)==NULL)
187 {
188 poly G=p_GcdMon(g,f,r);
189 if (!n_IsOne(pGetCoeff(G),r->cf)
190 || (!p_IsConstant(G,r)))
191 {
192 f=p_Div_mm(f,G,r);
193 g=p_Div_mm(g,G,r);
194 }
195 return G;
196 }
197 else if (pNext(f)==NULL)
198 {
199 poly G=p_GcdMon(f,g,r);
200 if (!n_IsOne(pGetCoeff(G),r->cf)
201 || (!p_IsConstant(G,r)))
202 {
203 f=p_Div_mm(f,G,r);
204 g=p_Div_mm(g,G,r);
205 }
206 return G;
207 }
208
209 Off(SW_RATIONAL);
210 CanonicalForm F,G,GCD;
211 if (rField_is_Q(r) || (rField_is_Zp(r))
212 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
213 {
214 bool b1=isOn(SW_USE_EZGCD_P);
215 setCharacteristic( rInternalChar(r) );
216 F=convSingPFactoryP( f,r );
217 G=convSingPFactoryP( g,r );
218 GCD=gcd(F,G);
219 if (!GCD.isOne())
220 {
221 p_Delete(&f,r);
222 p_Delete(&g,r);
223 if (getCharacteristic() == 0)
224 On (SW_RATIONAL);
225 F /= GCD;
226 G /= GCD;
227 if (getCharacteristic() == 0)
228 {
229 CanonicalForm denF= bCommonDen (F);
230 CanonicalForm denG= bCommonDen (G);
231 G *= denG;
232 F *= denF;
233 Off (SW_RATIONAL);
234 CanonicalForm gcddenFdenG= gcd (denG, denF);
235 denG /= gcddenFdenG;
236 denF /= gcddenFdenG;
237 On (SW_RATIONAL);
238 G *= denF;
239 F *= denG;
240 }
241 f=convFactoryPSingP( F, r);
242 g=convFactoryPSingP( G, r);
243 }
244 res=convFactoryPSingP( GCD , r);
245 if (!b1) Off (SW_USE_EZGCD_P);
246 }
247 // and over Q(a) / Fp(a)
248 else if ( r->cf->extRing )
249 {
250 if ( rField_is_Q_a(r)) setCharacteristic( 0 );
251 else setCharacteristic( rInternalChar(r) );
252 if (r->cf->extRing->qideal!=NULL)
253 {
254 bool b1=isOn(SW_USE_QGCD);
255 if ( rField_is_Q_a(r) ) On(SW_USE_QGCD);
256 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
257 r->cf->extRing);
258 Variable a=rootOf(mipo);
259 F=( convSingAPFactoryAP( f,a,r ) );
260 G=( convSingAPFactoryAP( g,a,r ) );
261 GCD=gcd(F,G);
262 if (!GCD.isOne())
263 {
264 p_Delete(&f,r);
265 p_Delete(&g,r);
266 if (getCharacteristic() == 0)
267 On (SW_RATIONAL);
268 F /= GCD;
269 G /= GCD;
270 if (getCharacteristic() == 0)
271 {
272 CanonicalForm denF= bCommonDen (F);
273 CanonicalForm denG= bCommonDen (G);
274 G *= denG;
275 F *= denF;
276 Off (SW_RATIONAL);
277 CanonicalForm gcddenFdenG= gcd (denG, denF);
278 denG /= gcddenFdenG;
279 denF /= gcddenFdenG;
280 On (SW_RATIONAL);
281 G *= denF;
282 F *= denG;
283 }
284 f= convFactoryAPSingAP( F,r );
285 g= convFactoryAPSingAP( G,r );
286 }
287 res= convFactoryAPSingAP( GCD,r );
288 prune (a);
289 if (!b1) Off(SW_USE_QGCD);
290 }
291 else
292 {
293 F=( convSingTrPFactoryP( f,r ) );
294 G=( convSingTrPFactoryP( g,r ) );
295 GCD=gcd(F,G);
296 if (!GCD.isOne())
297 {
298 p_Delete(&f,r);
299 p_Delete(&g,r);
300 if (getCharacteristic() == 0)
301 On (SW_RATIONAL);
302 F /= GCD;
303 G /= GCD;
304 if (getCharacteristic() == 0)
305 {
306 CanonicalForm denF= bCommonDen (F);
307 CanonicalForm denG= bCommonDen (G);
308 G *= denG;
309 F *= denF;
310 Off (SW_RATIONAL);
311 CanonicalForm gcddenFdenG= gcd (denG, denF);
312 denG /= gcddenFdenG;
313 denF /= gcddenFdenG;
314 On (SW_RATIONAL);
315 G *= denF;
316 F *= denG;
317 }
318 f= convFactoryPSingTrP( F,r );
319 g= convFactoryPSingTrP( G,r );
320 }
321 res= convFactoryPSingTrP( GCD,r );
322 }
323 }
324 else
325 WerrorS( feNotImplemented );
326 Off(SW_RATIONAL);
327 return res;
328}
329
330/*2 find the maximal exponent of var(i) in poly p*/
331int pGetExp_Var(poly p, int i, const ring r)
332{
333 int m=0;
334 int mm;
335 while (p!=NULL)
336 {
337 mm=p_GetExp(p,i,r);
338 if (mm>m) m=mm;
339 pIter(p);
340 }
341 return m;
342}
343
344// destroys f,g,x
345poly singclap_resultant ( poly f, poly g , poly x, const ring r)
346{
347 poly res=NULL;
348 int i=p_IsPurePower(x, r);
349 if (i==0)
350 {
351 WerrorS("3rd argument must be a ring variable");
352 goto resultant_returns_res;
353 }
354 if ((f==NULL) || (g==NULL))
355 goto resultant_returns_res;
356 // for now there is only the possibility to handle polynomials over
357 // Q and Fp ...
358 if (rField_is_Zp(r) || rField_is_Q(r)
359 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
360 {
361 Variable X(i);
362 setCharacteristic( rInternalChar(r) );
363 CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g,r ) );
364 res=convFactoryPSingP( resultant( F, G, X),r );
365 Off(SW_RATIONAL);
366 goto resultant_returns_res;
367 }
368 // and over Q(a) / Fp(a)
369 else if (r->cf->extRing!=NULL)
370 {
371 if (rField_is_Q_a(r)) setCharacteristic( 0 );
372 else setCharacteristic( rInternalChar(r) );
373 Variable X(i+rPar(r));
374 if (r->cf->extRing->qideal!=NULL)
375 {
376 //Variable X(i);
377 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
378 r->cf->extRing);
379 Variable a=rootOf(mipo);
380 CanonicalForm F( convSingAPFactoryAP( f,a,r ) ),
381 G( convSingAPFactoryAP( g,a,r ) );
382 res= convFactoryAPSingAP( resultant( F, G, X ),r );
383 prune (a);
384 }
385 else
386 {
387 //Variable X(i+rPar(currRing));
388 number nf,ng;
389 p_Cleardenom_n(f,r,nf);p_Cleardenom_n(g,r,ng);
390 int ef,eg;
391 ef=pGetExp_Var(f,i,r);
392 eg=pGetExp_Var(g,i,r);
393 CanonicalForm F( convSingTrPFactoryP( f,r ) ), G( convSingTrPFactoryP( g,r ) );
394 res= convFactoryPSingTrP( resultant( F, G, X ),r );
395 if ((nf!=NULL)&&(!n_IsOne(nf,r->cf)))
396 {
397 number n=n_Invers(nf,r->cf);
398 while(eg>0)
399 {
400 res=__p_Mult_nn(res,n,r);
401 eg--;
402 }
403 n_Delete(&n,r->cf);
404 }
405 n_Delete(&nf,r->cf);
406 if ((ng!=NULL)&&(!n_IsOne(ng,r->cf)))
407 {
408 number n=n_Invers(ng,r->cf);
409 while(ef>0)
410 {
411 res=__p_Mult_nn(res,n,r);
412 ef--;
413 }
414 n_Delete(&n,r->cf);
415 }
416 n_Delete(&ng,r->cf);
417 }
418 Off(SW_RATIONAL);
419 goto resultant_returns_res;
420 }
421 else
422 WerrorS( feNotImplemented );
423resultant_returns_res:
424 p_Delete(&f,r);
425 p_Delete(&g,r);
426 p_Delete(&x,r);
427 return res;
428}
429//poly singclap_resultant ( poly f, poly g , poly x)
430//{
431// int i=pVar(x);
432// if (i==0)
433// {
434// WerrorS("ringvar expected");
435// return NULL;
436// }
437// ideal I=idInit(1,1);
438//
439// // get the coeffs von f wrt. x:
440// I->m[0]=pCopy(f);
441// matrix ffi=mpCoeffs(I,i);
442// ffi->rank=1;
443// ffi->ncols=ffi->nrows;
444// ffi->nrows=1;
445// ideal fi=(ideal)ffi;
446//
447// // get the coeffs von g wrt. x:
448// I->m[0]=pCopy(g);
449// matrix ggi=mpCoeffs(I,i);
450// ggi->rank=1;
451// ggi->ncols=ggi->nrows;
452// ggi->nrows=1;
453// ideal gi=(ideal)ggi;
454//
455// // contruct the matrix:
456// int fn=IDELEMS(fi); //= deg(f,x)+1
457// int gn=IDELEMS(gi); //= deg(g,x)+1
458// matrix m=mpNew(fn+gn-2,fn+gn-2);
459// if(m==NULL)
460// {
461// return NULL;
462// }
463//
464// // enter the coeffs into m:
465// int j;
466// for(i=0;i<gn-1;i++)
467// {
468// for(j=0;j<fn;j++)
469// {
470// MATELEM(m,i+1,fn-j+i)=pCopy(fi->m[j]);
471// }
472// }
473// for(i=0;i<fn-1;i++)
474// {
475// for(j=0;j<gn;j++)
476// {
477// MATELEM(m,gn+i,gn-j+i)=pCopy(gi->m[j]);
478// }
479// }
480//
481// poly r=mpDet(m);
482//
483// idDelete(&fi);
484// idDelete(&gi);
485// idDelete((ideal *)&m);
486// return r;
487//}
488
489BOOLEAN singclap_extgcd ( poly f, poly g, poly &res, poly &pa, poly &pb , const ring r)
490{
491 // for now there is only the possibility to handle univariate
492 // polynomials over
493 // Q and Fp ...
494 res=NULL;pa=NULL;pb=NULL;
495 On(SW_SYMMETRIC_FF);
496 if ( rField_is_Q(r) || rField_is_Zp(r)
497 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
498 {
499 setCharacteristic( rInternalChar(r) );
500 CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g,r) );
501 CanonicalForm FpG=F+G;
502 if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
503 //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
504 {
505 Off(SW_RATIONAL);
506 WerrorS("not univariate");
507 return TRUE;
508 }
509 CanonicalForm Fa,Gb;
510 On(SW_RATIONAL);
511 res=convFactoryPSingP( extgcd( F, G, Fa, Gb ),r );
512 pa=convFactoryPSingP(Fa,r);
513 pb=convFactoryPSingP(Gb,r);
514 Off(SW_RATIONAL);
515 }
516 // and over Q(a) / Fp(a)
517 else if ( r->cf->extRing!=NULL )
518 {
519 if (rField_is_Q_a(r)) setCharacteristic( 0 );
520 else setCharacteristic( rInternalChar(r) );
521 CanonicalForm Fa,Gb;
522 if (r->cf->extRing->qideal!=NULL)
523 {
524 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
525 r->cf->extRing);
526 Variable a=rootOf(mipo);
527 CanonicalForm F( convSingAPFactoryAP( f,a,r ) ),
528 G( convSingAPFactoryAP( g,a,r ) );
529 CanonicalForm FpG=F+G;
530 if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
531 //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
532 {
533 WerrorS("not univariate");
534 return TRUE;
535 }
536 res= convFactoryAPSingAP( extgcd( F, G, Fa, Gb ),r );
537 pa=convFactoryAPSingAP(Fa,r);
538 pb=convFactoryAPSingAP(Gb,r);
539 prune (a);
540 }
541 else
542 {
543 CanonicalForm F( convSingTrPFactoryP( f, r ) ), G( convSingTrPFactoryP( g, r ) );
544 CanonicalForm FpG=F+G;
545 if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
546 //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
547 {
548 Off(SW_RATIONAL);
549 WerrorS("not univariate");
550 return TRUE;
551 }
552 res= convFactoryPSingTrP( extgcd( F, G, Fa, Gb ), r );
553 pa=convFactoryPSingTrP(Fa, r);
554 pb=convFactoryPSingTrP(Gb, r);
555 }
556 Off(SW_RATIONAL);
557 }
558 else
559 {
560 WerrorS( feNotImplemented );
561 return TRUE;
562 }
563#ifndef SING_NDEBUG
564 // checking the result of extgcd:
565 poly dummy;
566 dummy=p_Sub(p_Add_q(pp_Mult_qq(f,pa,r),pp_Mult_qq(g,pb,r),r),p_Copy(res,r),r);
567 if (dummy!=NULL)
568 {
569 PrintS("extgcd( ");p_Write(f,r);p_Write0(g,r);PrintS(" )\n");
570 PrintS("extgcd( ");p_Write(f,r);p_Write0(g,r);PrintS(" )\n");
571 p_Delete(&dummy,r);
572 }
573#endif
574 return FALSE;
575}
576
577poly singclap_pmult ( poly f, poly g, const ring r )
578{
579 poly res=NULL;
580 On(SW_RATIONAL);
581 if (rField_is_Zp(r) || rField_is_Q(r) || rField_is_Z(r)
582 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
583 {
584 if (rField_is_Z(r)) Off(SW_RATIONAL);
585 setCharacteristic( rInternalChar(r) );
586 CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g,r ) );
587 res = convFactoryPSingP( F * G,r );
588 }
589 else if (r->cf->extRing!=NULL)
590 {
591 if (rField_is_Q_a(r)) setCharacteristic( 0 );
592 else setCharacteristic( rInternalChar(r) );
593 if (r->cf->extRing->qideal!=NULL)
594 {
595 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
596 r->cf->extRing);
597 Variable a=rootOf(mipo);
598 CanonicalForm F( convSingAPFactoryAP( f,a,r ) ),
599 G( convSingAPFactoryAP( g,a,r ) );
600 res= convFactoryAPSingAP( F * G, r );
601 prune (a);
602 }
603 else
604 {
605 CanonicalForm F( convSingTrPFactoryP( f,r ) ), G( convSingTrPFactoryP( g,r ) );
606 res= convFactoryPSingTrP( F * G,r );
607 }
608 }
609#if 0 // not yet working
610 else if (rField_is_GF())
611 {
612 //Print("GF(%d^%d)\n",nfCharP,nfMinPoly[0]);
613 setCharacteristic( nfCharP,nfMinPoly[0], currRing->parameter[0][0] );
614 CanonicalForm F( convSingGFFactoryGF( f ) ), G( convSingGFFactoryGF( g ) );
615 res = convFactoryGFSingGF( F * G );
616 }
617#endif
618 else
619 WerrorS( feNotImplemented );
620 Off(SW_RATIONAL);
621 return res;
622}
623
624poly singclap_pdivide ( poly f, poly g, const ring r )
625{
626 poly res=NULL;
627
628 #ifdef HAVE_FLINT
629 #if __FLINT_RELEASE >= 20503
630 /*
631 If the division is not exact, control will pass to factory where the
632 polynomials can be divided using the ordering that factory chooses.
633 */
634 if (rField_is_Zp(r))
635 {
636 nmod_mpoly_ctx_t ctx;
637 if (!convSingRFlintR(ctx,r))
638 {
639 res = Flint_Divide_MP(f,0,g,0,ctx,r);
640 if (res != NULL)
641 return res;
642 }
643 }
644 else
645 if (rField_is_Q(r))
646 {
647 fmpq_mpoly_ctx_t ctx;
648 if (!convSingRFlintR(ctx,r))
649 {
650 res = Flint_Divide_MP(f,0,g,0,ctx,r);
651 if (res != NULL)
652 return res;
653 }
654 }
655 #endif
656 #endif
657
658 On(SW_RATIONAL);
659 if (rField_is_Zp(r) || rField_is_Q(r)
660 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
661 {
662 setCharacteristic( rInternalChar(r) );
663 CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g,r ) );
664 res = convFactoryPSingP( F / G,r );
665 }
666 // div is not implemented for ZZ coeffs in factory
667 else if (r->cf->extRing!=NULL)
668 {
669 if (rField_is_Q_a(r)) setCharacteristic( 0 );
670 else setCharacteristic( rInternalChar(r) );
671 if (r->cf->extRing->qideal!=NULL)
672 {
673 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
674 r->cf->extRing);
675 Variable a=rootOf(mipo);
676 CanonicalForm F( convSingAPFactoryAP( f,a,r ) ),
677 G( convSingAPFactoryAP( g,a,r ) );
678 res= convFactoryAPSingAP( F / G, r );
679 prune (a);
680 }
681 else
682 {
683 CanonicalForm F( convSingTrPFactoryP( f,r ) ), G( convSingTrPFactoryP( g,r ) );
684 res= convFactoryPSingTrP( F / G,r );
685 }
686 }
687#if 0 // not yet working
688 else if (rField_is_GF())
689 {
690 //Print("GF(%d^%d)\n",nfCharP,nfMinPoly[0]);
691 setCharacteristic( nfCharP,nfMinPoly[0], currRing->parameter[0][0] );
692 CanonicalForm F( convSingGFFactoryGF( f ) ), G( convSingGFFactoryGF( g ) );
693 res = convFactoryGFSingGF( F / G );
694 }
695#endif
696 else
697 WerrorS( feNotImplemented );
698 Off(SW_RATIONAL);
699 return res;
700}
701
702poly singclap_pmod ( poly f, poly g, const ring r )
703{
704 poly res=NULL;
705 On(SW_RATIONAL);
706 if (rField_is_Zp(r) || rField_is_Q(r)
707 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
708 {
709 setCharacteristic( rInternalChar(r) );
710 CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g,r ) );
711 CanonicalForm Q,R;
712 divrem(F,G,Q,R);
713 res = convFactoryPSingP(R,r);
714 //res = convFactoryPSingP( F-(F/G)*G,r );
715 }
716 // mod is not implemented for ZZ coeffs in factory
717 else if (r->cf->extRing!=NULL)
718 {
719 if (rField_is_Q_a(r)) setCharacteristic( 0 );
720 else setCharacteristic( rInternalChar(r) );
721 if (r->cf->extRing->qideal!=NULL)
722 {
723 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
724 r->cf->extRing);
725 Variable a=rootOf(mipo);
726 CanonicalForm F( convSingAPFactoryAP( f,a,r ) ),
727 G( convSingAPFactoryAP( g,a,r ) );
728 CanonicalForm Q,R;
729 divrem(F,G,Q,R);
730 res = convFactoryAPSingAP(R,r);
731 //res= convFactoryAPSingAP( F-(F/G)*G, r );
732 prune (a);
733 }
734 else
735 {
736 CanonicalForm F( convSingTrPFactoryP( f,r ) ), G( convSingTrPFactoryP( g,r ) );
737 CanonicalForm Q,R;
738 divrem(F,G,Q,R);
739 res = convFactoryPSingTrP(R,r);
740 //res= convFactoryPSingTrP( F-(F/G)*G,r );
741 }
742 }
743#if 0 // not yet working
744 else if (rField_is_GF())
745 {
746 //Print("GF(%d^%d)\n",nfCharP,nfMinPoly[0]);
747 setCharacteristic( nfCharP,nfMinPoly[0], currRing->parameter[0][0] );
748 CanonicalForm F( convSingGFFactoryGF( f ) ), G( convSingGFFactoryGF( g ) );
749 res = convFactoryGFSingGF( F / G );
750 }
751#endif
752 else
753 WerrorS( feNotImplemented );
754 Off(SW_RATIONAL);
755 return res;
756}
757
758#if 0
759// unused
760void singclap_divide_content ( poly f, const ring r )
761{
762 if ( f==NULL )
763 {
764 return;
765 }
766 else if ( pNext( f ) == NULL )
767 {
768 p_SetCoeff( f, n_Init( 1, r->cf ), r );
769 return;
770 }
771 else
772 {
773 if ( rField_is_Q_a(r) )
774 setCharacteristic( 0 );
775 else if ( rField_is_Zp_a(r) )
776 setCharacteristic( -rChar(r) );
777 else
778 return; /* not implemented*/
779
780 CFList L;
781 CanonicalForm g, h;
782 poly p = pNext(f);
783
784 // first attemp: find 2 smallest g:
785
786 number g1=pGetCoeff(f);
787 number g2=pGetCoeff(p); // p==pNext(f);
788 pIter(p);
789 int sz1=n_Size(g1, r->cf);
790 int sz2=n_Size(g2, r->cf);
791 if (sz1>sz2)
792 {
793 number gg=g1;
794 g1=g2; g2=gg;
795 int sz=sz1;
796 sz1=sz2; sz2=sz;
797 }
798 while (p!=NULL)
799 {
800 int n_sz=n_Size(pGetCoeff(p),r->cf);
801 if (n_sz<sz1)
802 {
803 sz2=sz1;
804 g2=g1;
805 g1=pGetCoeff(p);
806 sz1=n_sz;
807 if (sz1<=3) break;
808 }
809 else if(n_sz<sz2)
810 {
811 sz2=n_sz;
812 g2=pGetCoeff(p);
813 sz2=n_sz;
814 }
815 pIter(p);
816 }
817 g = convSingPFactoryP( NUM(((fraction)g1)), r->cf->extRing );
818 g = gcd( g, convSingPFactoryP( NUM(((fraction)g2)) , r->cf->extRing));
819
820 // second run: gcd's
821
822 p = f;
823 while ( (p != NULL) && (g != 1) && ( g != 0))
824 {
825 h = convSingPFactoryP( NUM(((fraction)pGetCoeff(p))), r->cf->extRing );
826 pIter( p );
827
828 g = gcd( g, h );
829
830 L.append( h );
831 }
832 if (( g == 1 ) || (g == 0))
833 {
834 // pTest(f);
835 return;
836 }
837 else
838 {
839 CFListIterator i;
840 for ( i = L, p = f; i.hasItem(); i++, p=pNext(p) )
841 {
842 fraction c=(fraction)pGetCoeff(p);
843 p_Delete(&(NUM(c)),r->cf->extRing); // 2nd arg used to be nacRing
844 NUM(c)=convFactoryPSingP( i.getItem() / g, r->cf->extRing );
845 //nTest((number)c);
846 //#ifdef LDEBUG
847 //number cn=(number)c;
848 //StringSetS(""); nWrite(nt); StringAppend(" ==> ");
849 //nWrite(cn);PrintS(StringEndS("\n")); // NOTE/TODO: use StringAppendS("\n"); omFree(s);
850 //#endif
851 }
852 }
853 // pTest(f);
854 }
855}
856#endif
857
858BOOLEAN count_Factors(ideal I, intvec *v,int j, poly &f, poly fac, const ring r)
859{
860 p_Test(f,r);
861 p_Test(fac,r);
862 int e=0;
863 if (!p_IsConstant(fac,r))
864 {
865#ifdef FACTORIZE2_DEBUG
866 printf("start count_Factors(%d), Fdeg=%d, factor deg=%d\n",j,p_Totaldegree(f,r),p_Totaldegree(fac,r));
867 p_wrp(fac,r);PrintLn();
868#endif
869 On(SW_RATIONAL);
870 CanonicalForm F, FAC,Q,R;
871 Variable a;
872 if (rField_is_Zp(r) || rField_is_Q(r)
873 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
874 {
875 F=convSingPFactoryP( f,r );
876 FAC=convSingPFactoryP( fac,r );
877 }
878 else if (r->cf->extRing!=NULL)
879 {
880 if (r->cf->extRing->qideal!=NULL)
881 {
882 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
883 r->cf->extRing);
884 a=rootOf(mipo);
885 F=convSingAPFactoryAP( f,a,r );
886 FAC=convSingAPFactoryAP( fac,a,r );
887 }
888 else
889 {
890 F=convSingTrPFactoryP( f,r );
891 FAC=convSingTrPFactoryP( fac,r );
892 }
893 }
894 else
895 WerrorS( feNotImplemented );
896
897 poly q;
898 loop
899 {
900 Q=F;
901 Q/=FAC;
902 R=Q;
903 R*=FAC;
904 R-=F;
905 if (R.isZero())
906 {
907 if (rField_is_Zp(r) || rField_is_Q(r)
908 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
909 {
910 q = convFactoryPSingP( Q,r );
911 }
912 else if (r->cf->extRing!=NULL)
913 {
914 if (r->cf->extRing->qideal!=NULL)
915 {
916 q= convFactoryAPSingAP( Q,r );
917 }
918 else
919 {
920 q= convFactoryPSingTrP( Q,r );
921 }
922 }
923 e++; p_Delete(&f,r); f=q; q=NULL; F=Q;
924 }
925 else
926 {
927 break;
928 }
929 }
930 if (r->cf->extRing!=NULL)
931 if (r->cf->extRing->qideal!=NULL)
932 prune (a);
933 if (e==0)
934 {
935 Off(SW_RATIONAL);
936 return FALSE;
937 }
938 }
939 else e=1;
940 I->m[j]=fac;
941 if (v!=NULL) (*v)[j]=e;
942 Off(SW_RATIONAL);
943 return TRUE;
944}
945
946VAR int singclap_factorize_retry;
947
948ideal singclap_factorize ( poly f, intvec ** v , int with_exps, const ring r)
949/* destroys f, sets *v */
950{
951 p_Test(f,r);
952#ifdef FACTORIZE2_DEBUG
953 printf("singclap_factorize, degree %ld\n",p_Totaldegree(f,r));
954#endif
955 // with_exps: 3,1 return only true factors, no exponents
956 // 2 return true factors and exponents
957 // 0 return coeff, factors and exponents
958 BOOLEAN save_errorreported=errorreported;
959
960 ideal res=NULL;
961
962 // handle factorize(0) =========================================
963 if (f==NULL)
964 {
965 res=idInit(1,1);
966 if (with_exps!=1)
967 {
968 (*v)=new intvec(1);
969 (**v)[0]=1;
970 }
971 return res;
972 }
973 // handle factorize(mon) =========================================
974 if (pNext(f)==NULL)
975 {
976 int i=0;
977 int n=0;
978 int e;
979 for(i=rVar(r);i>0;i--) if(p_GetExp(f,i,r)!=0) n++;
980 if (with_exps==0) n++; // with coeff
981 res=idInit(si_max(n,1),1);
982 switch(with_exps)
983 {
984 case 0: // with coef & exp.
985 res->m[0]=p_NSet(n_Copy(pGetCoeff(f),r->cf),r);
986 // no break
987 case 2: // with exp.
988 (*v)=new intvec(si_max(1,n));
989 (**v)[0]=1;
990 // no break
991 case 1: ;
992#ifdef TEST
993 default: ;
994#endif
995 }
996 if (n==0)
997 {
998 if (res->m[0]==NULL) res->m[0]=p_One(r);
999 // (**v)[0]=1; is already done
1000 }
1001 else
1002 {
1003 for(i=rVar(r);i>0;i--)
1004 {
1005 e=p_GetExp(f,i,r);
1006 if(e!=0)
1007 {
1008 n--;
1009 poly p=p_One(r);
1010 p_SetExp(p,i,1,r);
1011 p_Setm(p,r);
1012 res->m[n]=p;
1013 if (with_exps!=1) (**v)[n]=e;
1014 }
1015 }
1016 }
1017 p_Delete(&f,r);
1018 return res;
1019 }
1020 //PrintS("S:");p_Write(f,r);PrintLn();
1021 // use factory/libfac in general ==============================
1022 Variable dummy(-1); prune(dummy); // remove all (tmp.) extensions
1023 Off(SW_RATIONAL);
1024 On(SW_SYMMETRIC_FF);
1025 CFFList L;
1026 number N=NULL;
1027 number NN=NULL;
1028 number old_lead_coeff=n_Copy(pGetCoeff(f), r->cf);
1029
1030 Variable a;
1031 if (r->cf->convSingNFactoryN!=ndConvSingNFactoryN)
1032 {
1033 if (rField_is_Q(r) || rField_is_Q_a(r) || rField_is_Z(r)) /* Q, Q(a), Z */
1034 {
1035 //if (f!=NULL) // already tested at start of routine
1036 {
1037 number n0=n_Copy(pGetCoeff(f),r->cf);
1038 if (with_exps==0)
1039 N=n_Copy(n0,r->cf);
1040 if (rField_is_Z(r)) p_Content(f, r);
1041 else p_Cleardenom(f, r);
1042 //after here f should not have a denominator!! and no content
1043 //PrintS("S:");p_Write(f,r);PrintLn();
1044 NN=n_Div(n0,pGetCoeff(f),r->cf);
1045 n_Delete(&n0,r->cf);
1046 if (with_exps==0)
1047 {
1048 n_Delete(&N,r->cf);
1049 N=n_Copy(NN,r->cf);
1050 }
1051 }
1052 }
1053 else if (rField_is_Zp_a(r))
1054 {
1055 //if (f!=NULL) // already tested at start of routine
1056 if (singclap_factorize_retry==0)
1057 {
1058 number n0=n_Copy(pGetCoeff(f),r->cf);
1059 if (with_exps==0)
1060 N=n_Copy(n0,r->cf);
1061 p_Norm(f,r);
1062 p_Cleardenom(f, r);
1063 NN=n_Div(n0,pGetCoeff(f),r->cf);
1064 n_Delete(&n0,r->cf);
1065 if (with_exps==0)
1066 {
1067 n_Delete(&N,r->cf);
1068 N=n_Copy(NN,r->cf);
1069 }
1070 }
1071 }
1072 if (rField_is_Q(r) || rField_is_Zp(r) || rField_is_Z(r) || rField_is_Zn(r))
1073 {
1074 setCharacteristic( rInternalChar(r) );
1075 if (errorreported) goto notImpl; // char too large
1076 CanonicalForm F( convSingPFactoryP( f,r ) );
1077 L = factorize( F );
1078 }
1079 // and over Q(a) / Fp(a)
1080 else if (r->cf->extRing!=NULL)
1081 {
1082 if (rField_is_Q_a (r)) setCharacteristic (0);
1083 else setCharacteristic( rInternalChar(r) );
1084 if (errorreported) goto notImpl; // char too large
1085 if (r->cf->extRing->qideal!=NULL) /*algebraic extension */
1086 {
1087 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
1088 r->cf->extRing);
1089 a=rootOf(mipo);
1090 CanonicalForm F( convSingAPFactoryAP( f, a, r ) );
1091 L = factorize( F, a );
1092 prune(a);
1093 }
1094 else /* rational functions */
1095 {
1096 CanonicalForm F( convSingTrPFactoryP( f,r ) );
1097 L = factorize( F );
1098 }
1099 }
1100 else
1101 {
1102 goto notImpl;
1103 }
1104 }
1105 else
1106 {
1107 goto notImpl;
1108 }
1109 if (errorreported)
1110 {
1111 errorreported=FALSE;
1112 }
1113 {
1114 poly ff=p_Copy(f,r); // a copy for the retry stuff
1115 // the first factor should be a constant
1116 if ( ! L.getFirst().factor().inCoeffDomain() )
1117 L.insert(CFFactor(1,1));
1118 // convert into ideal
1119 int n = L.length();
1120 if (n==0) n=1;
1121 CFFListIterator J=L;
1122 int j=0;
1123 if (with_exps!=1)
1124 {
1125 if ((with_exps==2)&&(n>1))
1126 {
1127 n--;
1128 J++;
1129 }
1130 *v = new intvec( n );
1131 }
1132 res = idInit( n ,1);
1133 for ( ; J.hasItem(); J++, j++ )
1134 {
1135 if (with_exps!=1) (**v)[j] = J.getItem().exp();
1136 if (rField_is_Zp(r) || rField_is_Q(r)|| rField_is_Z(r)
1137 || rField_is_Zn(r)) /* Q, Fp, Z */
1138 {
1139 //count_Factors(res,*v,f, j, convFactoryPSingP( J.getItem().factor() );
1140 res->m[j] = convFactoryPSingP( J.getItem().factor(),r );
1141 }
1142#if 0
1143 else if (rField_is_GF())
1144 res->m[j] = convFactoryGFSingGF( J.getItem().factor() );
1145#endif
1146 else if (r->cf->extRing!=NULL) /* Q(a), Fp(a) */
1147 {
1148#ifndef SING_NDEBUG
1149 intvec *w=NULL;
1150 if (v!=NULL) w=*v;
1151#endif
1152 if (r->cf->extRing->qideal==NULL)
1153 {
1154#ifdef SING_NDEBUG
1155 res->m[j]= convFactoryPSingTrP( J.getItem().factor(),r );
1156#else
1157 if(!count_Factors(res,w,j,ff,convFactoryPSingTrP( J.getItem().factor(),r ),r))
1158 {
1159 if (w!=NULL)
1160 (*w)[j]=1;
1161 res->m[j]=p_One(r);
1162 }
1163#endif
1164 }
1165 else
1166 {
1167#ifdef SING_NDEBUG
1168 res->m[j]= convFactoryAPSingAP( J.getItem().factor(),r );
1169#else
1170 if (!count_Factors(res,w,j,ff,convFactoryAPSingAP( J.getItem().factor(),r ),r))
1171 {
1172 if (w!=NULL)
1173 (*w)[j]=1;
1174 res->m[j]=p_One(r);
1175 }
1176#endif
1177 }
1178 }
1179 }
1180 if (r->cf->extRing!=NULL)
1181 if (r->cf->extRing->qideal!=NULL)
1182 prune (a);
1183#ifndef SING_NDEBUG
1184 if ((r->cf->extRing!=NULL) && (!p_IsConstant(ff,r)))
1185 {
1186 singclap_factorize_retry++;
1187 if (singclap_factorize_retry<3)
1188 {
1189 int jj;
1190#ifdef FACTORIZE2_DEBUG
1191 printf("factorize_retry\n");
1192#endif
1193 intvec *ww=NULL;
1194 id_Test(res,r);
1195 ideal h=singclap_factorize ( ff, &ww , with_exps, r );
1196 id_Test(h,r);
1197 int l=(*v)->length();
1198 (*v)->resize(l+ww->length());
1199 for(jj=0;jj<ww->length();jj++)
1200 (**v)[jj+l]=(*ww)[jj];
1201 delete ww;
1202 ideal hh=idInit(IDELEMS(res)+IDELEMS(h),1);
1203 for(jj=IDELEMS(res)-1;jj>=0;jj--)
1204 {
1205 hh->m[jj]=res->m[jj];
1206 res->m[jj]=NULL;
1207 }
1208 for(jj=IDELEMS(h)-1;jj>=0;jj--)
1209 {
1210 hh->m[jj+IDELEMS(res)]=h->m[jj];
1211 h->m[jj]=NULL;
1212 }
1213 id_Delete(&res,r);
1214 id_Delete(&h,r);
1215 res=hh;
1216 id_Test(res,r);
1217 ff=NULL;
1218 }
1219 else
1220 {
1221 WarnS("problem with factorize");
1222#if 0
1223 pWrite(ff);
1224 idShow(res);
1225#endif
1226 id_Delete(&res,r);
1227 res=idInit(2,1);
1228 res->m[0]=p_One(r);
1229 res->m[1]=ff; ff=NULL;
1230 }
1231 }
1232#endif
1233 p_Delete(&ff,r);
1234 if (N!=NULL)
1235 {
1236 __p_Mult_nn(res->m[0],N,r);
1237 n_Delete(&N,r->cf);
1238 N=NULL;
1239 }
1240 // delete constants
1241 if (res!=NULL)
1242 {
1243 int i=IDELEMS(res)-1;
1244 int j=0;
1245 for(;i>=0;i--)
1246 {
1247 if ((res->m[i]!=NULL)
1248 && (pNext(res->m[i])==NULL)
1249 && (p_IsConstant(res->m[i],r)))
1250 {
1251 if (with_exps!=0)
1252 {
1253 p_Delete(&(res->m[i]),r);
1254 if ((v!=NULL) && ((*v)!=NULL))
1255 (**v)[i]=0;
1256 j++;
1257 }
1258 else if (i!=0)
1259 {
1260 while ((v!=NULL) && ((*v)!=NULL) && ((**v)[i]>1))
1261 {
1262 res->m[0]=p_Mult_q(res->m[0],p_Copy(res->m[i],r),r);
1263 (**v)[i]--;
1264 }
1265 res->m[0]=p_Mult_q(res->m[0],res->m[i],r);
1266 res->m[i]=NULL;
1267 if ((v!=NULL) && ((*v)!=NULL))
1268 (**v)[i]=1;
1269 j++;
1270 }
1271 }
1272 }
1273 if (j>0)
1274 {
1275 idSkipZeroes(res);
1276 if ((v!=NULL) && ((*v)!=NULL))
1277 {
1278 intvec *w=*v;
1279 int len=IDELEMS(res);
1280 *v = new intvec( len );
1281 for (i=0,j=0;i<si_min(w->length(),len);i++)
1282 {
1283 if((*w)[i]!=0)
1284 {
1285 (**v)[j]=(*w)[i]; j++;
1286 }
1287 }
1288 delete w;
1289 }
1290 }
1291 if (res->m[0]==NULL)
1292 {
1293 res->m[0]=p_One(r);
1294 }
1295 }
1296 }
1297 if (rField_is_Q_a(r) && (r->cf->extRing->qideal!=NULL))
1298 {
1299 int i=IDELEMS(res)-1;
1300 int stop=1;
1301 if (with_exps!=0) stop=0;
1302 for(;i>=stop;i--)
1303 {
1304 p_Norm(res->m[i],r);
1305 }
1306 if (with_exps==0) p_SetCoeff(res->m[0],old_lead_coeff,r);
1307 else n_Delete(&old_lead_coeff,r->cf);
1308 }
1309 else
1310 n_Delete(&old_lead_coeff,r->cf);
1311 errorreported=save_errorreported;
1312notImpl:
1313 prune(a);
1314 if (res==NULL)
1315 {
1316 WerrorS( feNotImplemented );
1317 if ((v!=NULL) && ((*v)!=NULL) &&(with_exps==2))
1318 {
1319 *v = new intvec( 1 );
1320 (*v)[0]=1;
1321 }
1322 res=idInit(2,1);
1323 res->m[0]=p_One(r);
1324 res->m[1]=f;
1325 }
1326 else p_Delete(&f,r);
1327 if (NN!=NULL)
1328 {
1329 n_Delete(&NN,r->cf);
1330 }
1331 if (N!=NULL)
1332 {
1333 n_Delete(&N,r->cf);
1334 }
1335 return res;
1336}
1337
1338ideal singclap_sqrfree ( poly f, intvec ** v , int with_exps, const ring r)
1339{
1340 p_Test(f,r);
1341#ifdef FACTORIZE2_DEBUG
1342 printf("singclap_sqrfree, degree %d\n",pTotaldegree(f));
1343#endif
1344 // with_exps: 3,1 return only true factors, no exponents
1345 // 2 return true factors and exponents
1346 // 0 return coeff, factors and exponents
1347 BOOLEAN save_errorreported=errorreported;
1348
1349 ideal res=NULL;
1350
1351 // handle factorize(0) =========================================
1352 if (f==NULL)
1353 {
1354 res=idInit(1,1);
1355 if (with_exps!=1 && with_exps!=3)
1356 {
1357 (*v)=new intvec(1);
1358 (**v)[0]=1;
1359 }
1360 return res;
1361 }
1362 // handle factorize(mon) =========================================
1363 if (pNext(f)==NULL)
1364 {
1365 int i=0;
1366 int n=0;
1367 int e;
1368 for(i=rVar(r);i>0;i--) if(p_GetExp(f,i,r)!=0) n++;
1369 if (with_exps==0 || with_exps==3) n++; // with coeff
1370 res=idInit(si_max(n,1),1);
1371 if(with_exps!=1)
1372 {
1373 (*v)=new intvec(si_max(1,n));
1374 (**v)[0]=1;
1375 }
1376 res->m[0]=p_NSet(n_Copy(pGetCoeff(f),r->cf),r);
1377 if (n==0)
1378 {
1379 res->m[0]=p_One(r);
1380 // (**v)[0]=1; is already done
1381 }
1382 else
1383 {
1384 for(i=rVar(r);i>0;i--)
1385 {
1386 e=p_GetExp(f,i,r);
1387 if(e!=0)
1388 {
1389 n--;
1390 poly p=p_One(r);
1391 p_SetExp(p,i,1,r);
1392 p_Setm(p,r);
1393 res->m[n]=p;
1394 if (with_exps!=1) (**v)[n]=e;
1395 }
1396 }
1397 }
1398 p_Delete(&f,r);
1399 return res;
1400 }
1401 //PrintS("S:");pWrite(f);PrintLn();
1402 // use factory/libfac in general ==============================
1403 Off(SW_RATIONAL);
1404 On(SW_SYMMETRIC_FF);
1405 CFFList L;
1406 number N=NULL;
1407 number NN=NULL;
1408 number old_lead_coeff=n_Copy(pGetCoeff(f), r->cf);
1409 Variable a;
1410
1411 if (!rField_is_Zp(r) && !rField_is_Zp_a(r)) /* Q, Q(a) */
1412 {
1413 //if (f!=NULL) // already tested at start of routine
1414 number n0=n_Copy(old_lead_coeff,r->cf);
1415 if (with_exps==0 || with_exps==3)
1416 N=n_Copy(n0,r->cf);
1417 p_Cleardenom(f, r);
1418 //after here f should not have a denominator!!
1419 //PrintS("S:");p_Write(f,r);PrintLn();
1420 NN=n_Div(n0,pGetCoeff(f),r->cf);
1421 n_Delete(&n0,r->cf);
1422 if (with_exps==0 || with_exps==3)
1423 {
1424 n_Delete(&N,r->cf);
1425 N=n_Copy(NN,r->cf);
1426 }
1427 }
1428 else if (rField_is_Zp_a(r))
1429 {
1430 //if (f!=NULL) // already tested at start of routine
1431 if (singclap_factorize_retry==0)
1432 {
1433 number n0=n_Copy(old_lead_coeff,r->cf);
1434 if (with_exps==0 || with_exps==3)
1435 N=n_Copy(n0,r->cf);
1436 p_Norm(f,r);
1437 p_Cleardenom(f, r);
1438 NN=n_Div(n0,pGetCoeff(f),r->cf);
1439 n_Delete(&n0,r->cf);
1440 if (with_exps==0 || with_exps==3)
1441 {
1442 n_Delete(&N,r->cf);
1443 N=n_Copy(NN,r->cf);
1444 }
1445 }
1446 }
1447 if (rField_is_Q(r) || rField_is_Zp(r)
1448 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
1449 {
1450 setCharacteristic( rInternalChar(r) );
1451 CanonicalForm F( convSingPFactoryP( f,r ) );
1452 L = sqrFree( F );
1453 }
1454 else if (r->cf->extRing!=NULL)
1455 {
1456 if (rField_is_Q_a (r)) setCharacteristic (0);
1457 else setCharacteristic( rInternalChar(r) );
1458 if (r->cf->extRing->qideal!=NULL)
1459 {
1460 CanonicalForm mipo=convSingPFactoryP(r->cf->extRing->qideal->m[0],
1461 r->cf->extRing);
1462 a=rootOf(mipo);
1463 CanonicalForm F( convSingAPFactoryAP( f, a, r ) );
1464 L= sqrFree (F);
1465 }
1466 else
1467 {
1468 CanonicalForm F( convSingTrPFactoryP( f,r ) );
1469 L = sqrFree( F );
1470 }
1471 }
1472 #if 0
1473 else if (rField_is_GF())
1474 {
1475 int c=rInternalChar(r);
1476 setCharacteristic( c, primepower(c) );
1477 CanonicalForm F( convSingGFFactoryGF( f ) );
1478 if (F.isUnivariate())
1479 {
1480 L = factorize( F );
1481 }
1482 else
1483 {
1484 goto notImpl;
1485 }
1486 }
1487 #endif
1488 else
1489 {
1490 goto notImpl;
1491 }
1492 {
1493 // convert into ideal
1494 int n = L.length();
1495 if (n==0) n=1;
1496 CFFListIterator J=L;
1497 int j=0;
1498 if (with_exps!=1)
1499 {
1500 if ((with_exps==2)&&(n>1))
1501 {
1502 n--;
1503 J++;
1504 }
1505 *v = new intvec( n );
1506 }
1507 else if (L.getFirst().factor().inCoeffDomain() && with_exps!=3)
1508 {
1509 n--;
1510 J++;
1511 }
1512 res = idInit( n ,1);
1513 for ( ; J.hasItem(); J++, j++ )
1514 {
1515 if (with_exps!=1 && with_exps!=3) (**v)[j] = J.getItem().exp();
1516 if (rField_is_Zp(r) || rField_is_Q(r)
1517 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
1518 res->m[j] = convFactoryPSingP( J.getItem().factor(),r );
1519 else if (r->cf->extRing!=NULL) /* Q(a), Fp(a) */
1520 {
1521 if (r->cf->extRing->qideal==NULL)
1522 res->m[j]=convFactoryPSingTrP( J.getItem().factor(),r );
1523 else
1524 res->m[j]=convFactoryAPSingAP( J.getItem().factor(),r );
1525 }
1526 }
1527 if (res->m[0]==NULL)
1528 {
1529 res->m[0]=p_One(r);
1530 }
1531 if (N!=NULL)
1532 {
1533 __p_Mult_nn(res->m[0],N,r);
1534 n_Delete(&N,r->cf);
1535 N=NULL;
1536 }
1537 }
1538 if (r->cf->extRing!=NULL)
1539 if (r->cf->extRing->qideal!=NULL)
1540 prune (a);
1541 if (rField_is_Q_a(r) && (r->cf->extRing->qideal!=NULL))
1542 {
1543 int i=IDELEMS(res)-1;
1544 int stop=1;
1545 if (with_exps!=0 || with_exps==3) stop=0;
1546 for(;i>=stop;i--)
1547 {
1548 p_Norm(res->m[i],r);
1549 }
1550 if (with_exps==0 || with_exps==3) p_SetCoeff(res->m[0],old_lead_coeff,r);
1551 else n_Delete(&old_lead_coeff,r->cf);
1552 }
1553 else
1554 n_Delete(&old_lead_coeff,r->cf);
1555 p_Delete(&f,r);
1556 errorreported=save_errorreported;
1557notImpl:
1558 if (res==NULL)
1559 WerrorS( feNotImplemented );
1560 if (NN!=NULL)
1561 {
1562 n_Delete(&NN,r->cf);
1563 }
1564 if (N!=NULL)
1565 {
1566 n_Delete(&N,r->cf);
1567 }
1568 return res;
1569}
1570
1571matrix singclap_irrCharSeries ( ideal I, const ring r)
1572{
1573 if (idIs0(I)) return mpNew(1,1);
1574
1575 // for now there is only the possibility to handle polynomials over
1576 // Q and Fp ...
1577 matrix res=NULL;
1578 int i;
1579 Off(SW_RATIONAL);
1580 On(SW_SYMMETRIC_FF);
1581 CFList L;
1582 ListCFList LL;
1583 if (rField_is_Q(r) || rField_is_Zp(r)
1584 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
1585 {
1586 setCharacteristic( rInternalChar(r) );
1587 for(i=0;i<IDELEMS(I);i++)
1588 {
1589 poly p=I->m[i];
1590 if (p!=NULL)
1591 {
1592 p=p_Copy(p,r);
1593 p_Cleardenom(p, r);
1594 L.append(convSingPFactoryP(p,r));
1595 p_Delete(&p,r);
1596 }
1597 }
1598 }
1599 // and over Q(a) / Fp(a)
1600 else if (nCoeff_is_transExt (r->cf))
1601 {
1602 setCharacteristic( rInternalChar(r) );
1603 for(i=0;i<IDELEMS(I);i++)
1604 {
1605 poly p=I->m[i];
1606 if (p!=NULL)
1607 {
1608 p=p_Copy(p,r);
1609 p_Cleardenom(p, r);
1610 L.append(convSingTrPFactoryP(p,r));
1611 p_Delete(&p,r);
1612 }
1613 }
1614 }
1615 else
1616 {
1617 WerrorS( feNotImplemented );
1618 return res;
1619 }
1620
1621 // a very bad work-around --- FIX IT in libfac
1622 // should be fixed as of 2001/6/27
1623 int tries=0;
1624 int m,n;
1625 ListIterator<CFList> LLi;
1626 loop
1627 {
1628 LL=irrCharSeries(L);
1629 m= LL.length(); // Anzahl Zeilen
1630 n=0;
1631 for ( LLi = LL; LLi.hasItem(); LLi++ )
1632 {
1633 n = si_max(LLi.getItem().length(),n);
1634 }
1635 if ((m!=0) && (n!=0)) break;
1636 tries++;
1637 if (tries>=5) break;
1638 }
1639 if ((m==0) || (n==0))
1640 {
1641 Warn("char_series returns %d x %d matrix from %d input polys (%d)",
1642 m,n,IDELEMS(I)+1,LL.length());
1643 iiWriteMatrix((matrix)I,"I",2,r,0);
1644 m=si_max(m,1);
1645 n=si_max(n,1);
1646 }
1647 res=mpNew(m,n);
1648 CFListIterator Li;
1649 for ( m=1, LLi = LL; LLi.hasItem(); LLi++, m++ )
1650 {
1651 for (n=1, Li = LLi.getItem(); Li.hasItem(); Li++, n++)
1652 {
1653 if (rField_is_Q(r) || rField_is_Zp(r)
1654 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
1655 MATELEM(res,m,n)=convFactoryPSingP(Li.getItem(),r);
1656 else
1657 MATELEM(res,m,n)=convFactoryPSingTrP(Li.getItem(),r);
1658 }
1659 }
1660 Off(SW_RATIONAL);
1661 return res;
1662}
1663
1664char* singclap_neworder ( ideal I, const ring r)
1665{
1666 int i;
1667 Off(SW_RATIONAL);
1668 On(SW_SYMMETRIC_FF);
1669 CFList L;
1670 if (rField_is_Q(r) || rField_is_Zp(r)
1671 || (rField_is_Zn(r)&&(r->cf->convSingNFactoryN!=ndConvSingNFactoryN)))
1672 {
1673 setCharacteristic( rInternalChar(r) );
1674 for(i=0;i<IDELEMS(I);i++)
1675 {
1676 poly p=I->m[i];
1677 if (p!=NULL)
1678 {
1679 p=p_Copy(p,r);
1680 p_Cleardenom(p, r);
1681 L.append(convSingPFactoryP(p,r));
1682 }
1683 }
1684 }
1685 // and over Q(a) / Fp(a)
1686 else if (nCoeff_is_transExt (r->cf))
1687 {
1688 setCharacteristic( rInternalChar(r) );
1689 for(i=0;i<IDELEMS(I);i++)
1690 {
1691 poly p=I->m[i];
1692 if (p!=NULL)
1693 {
1694 p=p_Copy(p,r);
1695 p_Cleardenom(p, r);
1696 L.append(convSingTrPFactoryP(p,r));
1697 }
1698 }
1699 }
1700 else
1701 {
1702 WerrorS( feNotImplemented );
1703 return NULL;
1704 }
1705
1706 List<int> IL=neworderint(L);
1707 ListIterator<int> Li;
1708 StringSetS("");
1709 Li = IL;
1710 int offs=rPar(r);
1711 int* mark=(int*)omAlloc0((rVar(r)+offs)*sizeof(int));
1712 int cnt=rVar(r)+offs;
1713 loop
1714 {
1715 if(! Li.hasItem()) break;
1716 BOOLEAN done=TRUE;
1717 i=Li.getItem()-1;
1718 mark[i]=1;
1719 if (i<offs)
1720 {
1721 done=FALSE;
1722 //StringAppendS(r->parameter[i]);
1723 }
1724 else
1725 {
1726 StringAppendS(r->names[i-offs]);
1727 }
1728 Li++;
1729 cnt--;
1730 if(cnt==0) break;
1731 if (done) StringAppendS(",");
1732 }
1733 for(i=0;i<rVar(r)+offs;i++)
1734 {
1735 BOOLEAN done=TRUE;
1736 if(mark[i]==0)
1737 {
1738 if (i<offs)
1739 {
1740 done=FALSE;
1741 //StringAppendS(r->parameter[i]);
1742 }
1743 else
1744 {
1745 StringAppendS(r->names[i-offs]);
1746 }
1747 cnt--;
1748 if(cnt==0) break;
1749 if (done) StringAppendS(",");
1750 }
1751 }
1752 char * s=StringEndS();
1753 if (s[strlen(s)-1]==',') s[strlen(s)-1]='\0';
1754 return s;
1755}
1756
1757poly singclap_det( const matrix m, const ring s )
1758{
1759 int r=m->rows();
1760 if (r!=m->cols())
1761 {
1762 Werror("det of %d x %d matrix",r,m->cols());
1763 return NULL;
1764 }
1765 poly res=NULL;
1766 CFMatrix M(r,r);
1767 int i,j;
1768 for(i=r;i>0;i--)
1769 {
1770 for(j=r;j>0;j--)
1771 {
1772 M(i,j)=convSingPFactoryP(MATELEM(m,i,j),s);
1773 }
1774 }
1775 res= convFactoryPSingP( determinant(M,r),s ) ;
1776 Off(SW_RATIONAL);
1777 return res;
1778}
1779
1780int singclap_det_i( intvec * m, const ring /*r*/)
1781{
1782// assume( r == currRing ); // Anything else is not guaranted to work!
1783
1784 setCharacteristic( 0 ); // ?
1785 CFMatrix M(m->rows(),m->cols());
1786 int i,j;
1787 for(i=m->rows();i>0;i--)
1788 {
1789 for(j=m->cols();j>0;j--)
1790 {
1791 M(i,j)=IMATELEM(*m,i,j);
1792 }
1793 }
1794 int res= convFactoryISingI( determinant(M,m->rows()) ) ;
1795 return res;
1796}
1797
1798number singclap_det_bi( bigintmat * m, const coeffs cf)
1799{
1800 assume(m->basecoeffs()==cf);
1801 CFMatrix M(m->rows(),m->cols());
1802 int i,j;
1803 BOOLEAN setchar=TRUE;
1804 for(i=m->rows();i>0;i--)
1805 {
1806 for(j=m->cols();j>0;j--)
1807 {
1808 M(i,j)=n_convSingNFactoryN(BIMATELEM(*m,i,j),setchar,cf);
1809 setchar=FALSE;
1810 }
1811 }
1812 number res=n_convFactoryNSingN( determinant(M,m->rows()),cf ) ;
1813 return res;
1814}
1815
1816#if defined(HAVE_NTL) || defined(AHVE_FLINT)
1817matrix singntl_HNF(matrix m, const ring s )
1818{
1819 int r=m->rows();
1820 if (r!=m->cols())
1821 {
1822 Werror("HNF of %d x %d matrix",r,m->cols());
1823 return NULL;
1824 }
1825
1826 matrix res=mpNew(r,r);
1827
1828 if (rField_is_Q(s))
1829 {
1830
1831 CFMatrix M(r,r);
1832 int i,j;
1833 for(i=r;i>0;i--)
1834 {
1835 for(j=r;j>0;j--)
1836 {
1837 M(i,j)=convSingPFactoryP(MATELEM(m,i,j),s );
1838 }
1839 }
1840 CFMatrix *MM=cf_HNF(M);
1841 for(i=r;i>0;i--)
1842 {
1843 for(j=r;j>0;j--)
1844 {
1845 MATELEM(res,i,j)=convFactoryPSingP((*MM)(i,j),s);
1846 }
1847 }
1848 delete MM;
1849 }
1850 return res;
1851}
1852
1853intvec* singntl_HNF(intvec* m)
1854{
1855 int r=m->rows();
1856 if (r!=m->cols())
1857 {
1858 Werror("HNF of %d x %d matrix",r,m->cols());
1859 return NULL;
1860 }
1861 setCharacteristic( 0 );
1862 CFMatrix M(r,r);
1863 int i,j;
1864 for(i=r;i>0;i--)
1865 {
1866 for(j=r;j>0;j--)
1867 {
1868 M(i,j)=IMATELEM(*m,i,j);
1869 }
1870 }
1871 CFMatrix *MM=cf_HNF(M);
1872 intvec *mm=ivCopy(m);
1873 for(i=r;i>0;i--)
1874 {
1875 for(j=r;j>0;j--)
1876 {
1877 IMATELEM(*mm,i,j)=convFactoryISingI((*MM)(i,j));
1878 }
1879 }
1880 delete MM;
1881 return mm;
1882}
1883
1884bigintmat* singntl_HNF(bigintmat* b)
1885{
1886 int r=b->rows();
1887 if (r!=b->cols())
1888 {
1889 Werror("HNF of %d x %d matrix",r,b->cols());
1890 return NULL;
1891 }
1892 setCharacteristic( 0 );
1893 CFMatrix M(r,r);
1894 int i,j;
1895 for(i=r;i>0;i--)
1896 {
1897 for(j=r;j>0;j--)
1898 {
1899 M(i,j)=n_convSingNFactoryN(BIMATELEM(*b,i,j),FALSE,b->basecoeffs());
1900 }
1901 }
1902 CFMatrix *MM=cf_HNF(M);
1903 bigintmat *mm=bimCopy(b);
1904 for(i=r;i>0;i--)
1905 {
1906 for(j=r;j>0;j--)
1907 {
1908 BIMATELEM(*mm,i,j)=n_convFactoryNSingN((*MM)(i,j),b->basecoeffs());
1909 }
1910 }
1911 delete MM;
1912 return mm;
1913}
1914
1915matrix singntl_LLL(matrix m, const ring s )
1916{
1917 int r=m->rows();
1918 int c=m->cols();
1919 matrix res=mpNew(r,c);
1920 if (rField_is_Q(s))
1921 {
1922 CFMatrix M(r,c);
1923 int i,j;
1924 for(i=r;i>0;i--)
1925 {
1926 for(j=c;j>0;j--)
1927 {
1928 M(i,j)=convSingPFactoryP(MATELEM(m,i,j),s);
1929 }
1930 }
1931 CFMatrix *MM=cf_LLL(M);
1932 for(i=r;i>0;i--)
1933 {
1934 for(j=c;j>0;j--)
1935 {
1936 MATELEM(res,i,j)=convFactoryPSingP((*MM)(i,j),s);
1937 }
1938 }
1939 delete MM;
1940 }
1941 return res;
1942}
1943
1944intvec* singntl_LLL(intvec* m)
1945{
1946 int r=m->rows();
1947 int c=m->cols();
1948 setCharacteristic( 0 );
1949 CFMatrix M(r,c);
1950 int i,j;
1951 for(i=r;i>0;i--)
1952 {
1953 for(j=c;j>0;j--)
1954 {
1955 M(i,j)=IMATELEM(*m,i,j);
1956 }
1957 }
1958 CFMatrix *MM=cf_LLL(M);
1959 intvec *mm=ivCopy(m);
1960 for(i=r;i>0;i--)
1961 {
1962 for(j=c;j>0;j--)
1963 {
1964 IMATELEM(*mm,i,j)=convFactoryISingI((*MM)(i,j));
1965 }
1966 }
1967 delete MM;
1968 return mm;
1969}
1970#else
1971matrix singntl_HNF(matrix m, const ring s )
1972{
1973 WerrorS("NTL/FLINT missing");
1974 return NULL;
1975}
1976
1977intvec* singntl_HNF(intvec* m)
1978{
1979 WerrorS("NTL/FLINT missing");
1980 return NULL;
1981}
1982
1983matrix singntl_LLL(matrix m, const ring s )
1984{
1985 WerrorS("NTL/FLINT missing");
1986 return NULL;
1987}
1988
1989intvec* singntl_LLL(intvec* m)
1990{
1991 WerrorS("NTL/FLINT missing");
1992 return NULL;
1993}
1994#endif
1995
1996#ifdef HAVE_NTL
1997matrix singntl_rref(matrix m, const ring R)
1998{
1999 int r=m->rows();
2000 int c=m->cols();
2001 int i,j;
2002 matrix M=mpNew(r,c);
2003 if (rField_is_Zp(R))
2004 {
2005 zz_p::init(rChar(R));
2006 mat_zz_p *NTLM=new mat_zz_p;
2007 NTLM->SetDims(r,c);
2008 for(i=r;i>0;i--)
2009 {
2010 for(j=c;j>0;j--)
2011 {
2012 poly h=MATELEM(m,i,j);
2013 if (h!=NULL)
2014 {
2015 if (p_Totaldegree(h,R)==0)
2016 {
2017 (*NTLM)(i,j)=(long)p_GetCoeff(h,R);
2018 }
2019 else
2020 {
2021 WerrorS("smatrix for rref is not constant");
2022 return M;
2023 }
2024 }
2025 }
2026 }
2027 gauss(*NTLM);
2028 for(i=r;i>0;i--)
2029 {
2030 for(j=c;j>0;j--)
2031 {
2032 number n=n_Init(rep((*NTLM)(i,j)),R->cf);
2033 if(!n_IsZero(n,R->cf))
2034 {
2035 poly p=p_NSet(n,R);
2036 MATELEM(M,i,j)=p;
2037 }
2038 }
2039 }
2040 delete NTLM;
2041 }
2042 else
2043 {
2044 WerrorS("not implemented for these coefficients");
2045 }
2046 return M;
2047}
2048#endif
2049
2050#ifdef HAVE_NTL
2051ideal singntl_rref(ideal m, const ring R) /*assume smatrix m*/
2052{
2053 int r=m->rank;
2054 int c=m->ncols;
2055 int i,j;
2056 ideal M=idInit(c,r);
2057 if (rField_is_Zp(R))
2058 {
2059 zz_p::init(rChar(R));
2060 mat_zz_p *NTLM=new mat_zz_p;
2061 NTLM->SetDims(r,c);
2062 for(j=c-1;j>=0;j--)
2063 {
2064 poly h=m->m[j];
2065 while(h!=NULL)
2066 {
2067 i=p_GetComp(h,R);
2068 if (p_Totaldegree(h,R)==0)
2069 (*NTLM)(i,j+1)=(long)p_GetCoeff(h,R);
2070 else
2071 {
2072 WerrorS("smatrix for rref is not constant");
2073 return M;
2074 }
2075 pIter(h);
2076 }
2077 }
2078 gauss(*NTLM);
2079 for(i=r;i>0;i--)
2080 {
2081 for(j=c;j>0;j--)
2082 {
2083 number n=n_Init(rep((*NTLM)(i,j)),R->cf);
2084 if(!n_IsZero(n,R->cf))
2085 {
2086 poly p=p_NSet(n,R);
2087 p_SetComp(p,i,R);
2088 M->m[j-1]=p_Add_q(M->m[j-1],p,R);
2089 }
2090 }
2091 }
2092 delete NTLM;
2093 }
2094 else
2095 {
2096 WerrorS("not implemented for these coefficients");
2097 }
2098 return M;
2099}
2100#endif
2101
2102#if defined(HAVE_NTL) || defined(HAVE_FLINT)
2103ideal singclap_absFactorize ( poly f, ideal & mipos, intvec ** exps, int & numFactors, const ring r)
2104{
2105 p_Test(f, r);
2106
2107 ideal res=NULL;
2108
2109 int offs = rPar(r);
2110 if (f==NULL)
2111 {
2112 res= idInit (1, 1);
2113 mipos= idInit (1, 1);
2114 mipos->m[0]= convFactoryPSingTrP (Variable (offs), r); //overkill
2115 (*exps)=new intvec (1);
2116 (**exps)[0]= 1;
2117 numFactors= 0;
2118 return res;
2119 }
2120 CanonicalForm F( convSingTrPFactoryP( f, r) );
2121
2122 bool isRat= isOn (SW_RATIONAL);
2123 if (!isRat)
2124 On (SW_RATIONAL);
2125
2126 CFAFList absFactors= absFactorize (F);
2127
2128 int n= absFactors.length();
2129 *exps=new intvec (n);
2130
2131 res= idInit (n, 1);
2132
2133 mipos= idInit (n, 1);
2134
2135 Variable x= Variable (offs);
2136 Variable alpha;
2137 int i= 0;
2138 numFactors= 0;
2139 int count;
2140 CFAFListIterator iter= absFactors;
2141 CanonicalForm lead= iter.getItem().factor();
2142 if (iter.getItem().factor().inCoeffDomain())
2143 {
2144 i++;
2145 iter++;
2146 }
2147 for (; iter.hasItem(); iter++, i++)
2148 {
2149 (**exps)[i]= iter.getItem().exp();
2150 alpha= iter.getItem().minpoly().mvar();
2151 if (iter.getItem().minpoly().isOne())
2152 lead /= power (bCommonDen (iter.getItem().factor()), iter.getItem().exp());
2153 else
2154 lead /= power (power (bCommonDen (iter.getItem().factor()), degree (iter.getItem().minpoly())), iter.getItem().exp());
2155 res->m[i]= convFactoryPSingTrP (replacevar (iter.getItem().factor()*bCommonDen (iter.getItem().factor()), alpha, x), r);
2156 if (iter.getItem().minpoly().isOne())
2157 {
2158 count= iter.getItem().exp();
2159 mipos->m[i]= convFactoryPSingTrP (x,r);
2160 }
2161 else
2162 {
2163 count= iter.getItem().exp()*degree (iter.getItem().minpoly());
2164 mipos->m[i]= convFactoryPSingTrP (replacevar (iter.getItem().minpoly(), alpha, x), r);
2165 }
2166 if (!iter.getItem().minpoly().isOne())
2167 prune (alpha);
2168 numFactors += count;
2169 }
2170 if (!isRat)
2171 Off (SW_RATIONAL);
2172
2173 (**exps)[0]= 1;
2174 res->m[0]= convFactoryPSingTrP (lead, r);
2175 mipos->m[0]= convFactoryPSingTrP (x, r);
2176 return res;
2177}
2178
2179#else
2180ideal singclap_absFactorize ( poly f, ideal & mipos, intvec ** exps, int & numFactors, const ring r)
2181{
2182 WerrorS("NTL/FLINT missing: absFactorize");
2183 return NULL;
2184}
2185
2186#endif /* HAVE_NTL */
2187
2188int * Zp_roots(poly p, const ring r)
2189{
2190 CanonicalForm pp=convSingPFactoryP(p,r);
2191 return Zp_roots(pp);
2192}
auxiliary.h
All the auxiliary stuff.
si_max
static int si_max(const int a, const int b)
Definition: auxiliary.h:124
BOOLEAN
int BOOLEAN
Definition: auxiliary.h:87
TRUE
#define TRUE
Definition: auxiliary.h:100
FALSE
#define FALSE
Definition: auxiliary.h:96
si_min
static int si_min(const int a, const int b)
Definition: auxiliary.h:125
bimCopy
bigintmat * bimCopy(const bigintmat *b)
same as copy constructor - apart from it being able to accept NULL as input
Definition: bigintmat.cc:405
BIMATELEM
#define BIMATELEM(M, I, J)
Definition: bigintmat.h:133
divrem
void divrem(const CanonicalForm &f, const CanonicalForm &g, CanonicalForm &q, CanonicalForm &r)
Definition: canonicalform.cc:1032
isOn
bool isOn(int sw)
switches
Definition: canonicalform.cc:1971
On
void On(int sw)
switches
Definition: canonicalform.cc:1957
Off
void Off(int sw)
switches
Definition: canonicalform.cc:1964
power
CanonicalForm power(const CanonicalForm &f, int n)
exponentiation
Definition: canonicalform.cc:1896
replacevar
CanonicalForm FACTORY_PUBLIC replacevar(const CanonicalForm &, const Variable &, const Variable &)
CanonicalForm replacevar ( const CanonicalForm & f, const Variable & x1, const Variable & x2 )
Definition: cf_ops.cc:271
degree
int degree(const CanonicalForm &f)
Definition: canonicalform.h:316
pp
CanonicalForm FACTORY_PUBLIC pp(const CanonicalForm &)
CanonicalForm pp ( const CanonicalForm & f )
Definition: cf_gcd.cc:676
setCharacteristic
void FACTORY_PUBLIC setCharacteristic(int c)
Definition: cf_char.cc:28
CFFactor
Factor< CanonicalForm > CFFactor
Definition: canonicalform.h:392
getCharacteristic
int FACTORY_PUBLIC getCharacteristic()
Definition: cf_char.cc:70
irrCharSeries
ListCFList irrCharSeries(const CFList &PS)
irreducible characteristic series
Definition: cfCharSets.cc:568
neworderint
IntList neworderint(const CFList &PolyList)
Definition: cfCharSets.cc:88
N
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:56
l
int l
Definition: cfEzgcd.cc:100
m
int m
Definition: cfEzgcd.cc:128
i
int i
Definition: cfEzgcd.cc:132
x
Variable x
Definition: cfModGcd.cc:4082
p
int p
Definition: cfModGcd.cc:4078
g
g
Definition: cfModGcd.cc:4090
cf
CanonicalForm cf
Definition: cfModGcd.cc:4083
b
CanonicalForm b
Definition: cfModGcd.cc:4103
extgcd
CanonicalForm extgcd(const CanonicalForm &f, const CanonicalForm &g, CanonicalForm &a, CanonicalForm &b)
CanonicalForm extgcd ( const CanonicalForm & f, const CanonicalForm & g, CanonicalForm & a,...
Definition: cfUnivarGcd.cc:174
bCommonDen
CanonicalForm bCommonDen(const CanonicalForm &f)
CanonicalForm bCommonDen ( const CanonicalForm & f )
Definition: cf_algorithm.cc:293
determinant
CanonicalForm FACTORY_PUBLIC determinant(const CFMatrix &M, int n)
Definition: cf_linsys.cc:222
sqrFree
CFFList FACTORY_PUBLIC sqrFree(const CanonicalForm &f, bool sort=false)
squarefree factorization
Definition: cf_factor.cc:957
factorize
CFFList FACTORY_PUBLIC factorize(const CanonicalForm &f, bool issqrfree=false)
factorization over or
Definition: cf_factor.cc:405
resultant
CanonicalForm FACTORY_PUBLIC resultant(const CanonicalForm &f, const CanonicalForm &g, const Variable &x)
CanonicalForm resultant ( const CanonicalForm & f, const CanonicalForm & g, const Variable & x )
Definition: cf_resultant.cc:173
SW_USE_QGCD
static const int SW_USE_QGCD
set to 1 to use Encarnacion GCD over Q(a)
Definition: cf_defs.h:43
SW_USE_EZGCD_P
static const int SW_USE_EZGCD_P
set to 1 to use EZGCD over F_q
Definition: cf_defs.h:37
SW_RATIONAL
static const int SW_RATIONAL
set to 1 for computations over Q
Definition: cf_defs.h:31
SW_SYMMETRIC_FF
static const int SW_SYMMETRIC_FF
set to 1 for symmetric representation over F_q
Definition: cf_defs.h:33
cf_LLL
CFMatrix * cf_LLL(CFMatrix &A)
performs LLL reduction.
Definition: cf_hnf.cc:66
cf_HNF
CFMatrix * cf_HNF(CFMatrix &A)
The input matrix A is an n x m matrix of rank m (so n >= m), and D is a multiple of the determinant o...
Definition: cf_hnf.cc:44
cf_roots.h
f
FILE * f
Definition: checklibs.c:9
convSingPFactoryP
CanonicalForm convSingPFactoryP(poly p, const ring r)
Definition: clapconv.cc:136
convSingTrP
BOOLEAN convSingTrP(poly p, const ring r)
Definition: clapconv.cc:352
convFactoryPSingTrP
poly convFactoryPSingTrP(const CanonicalForm &f, const ring r)
Definition: clapconv.cc:366
convFactoryAPSingAP
poly convFactoryAPSingAP(const CanonicalForm &f, const ring r)
Definition: clapconv.cc:184
convSingAPFactoryAP
CanonicalForm convSingAPFactoryAP(poly p, const Variable &a, const ring r)
Definition: clapconv.cc:148
convFactoryPSingP
poly convFactoryPSingP(const CanonicalForm &f, const ring r)
Definition: clapconv.cc:40
convFactoryISingI
int convFactoryISingI(const CanonicalForm &f)
Definition: clapconv.cc:142
convSingTrPFactoryP
CanonicalForm convSingTrPFactoryP(poly p, const ring r)
Definition: clapconv.cc:316
singclap_resultant
poly singclap_resultant(poly f, poly g, poly x, const ring r)
Definition: clapsing.cc:345
count_Factors
BOOLEAN count_Factors(ideal I, intvec *v, int j, poly &f, poly fac, const ring r)
Definition: clapsing.cc:858
singntl_HNF
matrix singntl_HNF(matrix m, const ring s)
Definition: clapsing.cc:1817
singclap_pmod
poly singclap_pmod(poly f, poly g, const ring r)
Definition: clapsing.cc:702
singclap_factorize
ideal singclap_factorize(poly f, intvec **v, int with_exps, const ring r)
Definition: clapsing.cc:948
singclap_irrCharSeries
matrix singclap_irrCharSeries(ideal I, const ring r)
Definition: clapsing.cc:1571
pGetExp_Var
int pGetExp_Var(poly p, int i, const ring r)
Definition: clapsing.cc:331
singclap_factorize_retry
VAR int singclap_factorize_retry
Definition: clapsing.cc:946
singclap_pdivide
poly singclap_pdivide(poly f, poly g, const ring r)
Definition: clapsing.cc:624
Zp_roots
int * Zp_roots(poly p, const ring r)
Definition: clapsing.cc:2188
out_cf
void out_cf(const char *s1, const CanonicalForm &f, const char *s2)
Definition: cf_factor.cc:99
singclap_extgcd
BOOLEAN singclap_extgcd(poly f, poly g, poly &res, poly &pa, poly &pb, const ring r)
Definition: clapsing.cc:489
singntl_rref
matrix singntl_rref(matrix m, const ring R)
Definition: clapsing.cc:1997
singclap_det
poly singclap_det(const matrix m, const ring s)
Definition: clapsing.cc:1757
singclap_gcd_and_divide
poly singclap_gcd_and_divide(poly &f, poly &g, const ring r)
clears denominators of f and g, divides by gcd(f,g)
Definition: clapsing.cc:170
singntl_LLL
matrix singntl_LLL(matrix m, const ring s)
Definition: clapsing.cc:1915
singclap_det_bi
number singclap_det_bi(bigintmat *m, const coeffs cf)
Definition: clapsing.cc:1798
singclap_pmult
poly singclap_pmult(poly f, poly g, const ring r)
Definition: clapsing.cc:577
singclap_gcd_r
poly singclap_gcd_r(poly f, poly g, const ring r)
Definition: clapsing.cc:68
singclap_det_i
int singclap_det_i(intvec *m, const ring)
Definition: clapsing.cc:1780
singclap_absFactorize
ideal singclap_absFactorize(poly f, ideal &mipos, intvec **exps, int &numFactors, const ring r)
Definition: clapsing.cc:2103
singclap_sqrfree
ideal singclap_sqrfree(poly f, intvec **v, int with_exps, const ring r)
Definition: clapsing.cc:1338
singclap_neworder
char * singclap_neworder(ideal I, const ring r)
Definition: clapsing.cc:1664
CanonicalForm
factory's main class
Definition: canonicalform.h:86
CanonicalForm::isOne
CF_NO_INLINE bool isOne() const
CanonicalForm::inCoeffDomain
bool inCoeffDomain() const
Definition: canonicalform.cc:122
CanonicalForm::isUnivariate
bool isUnivariate() const
Definition: canonicalform.cc:155
ListIterator::getItem
T & getItem() const
Definition: ftmpl_list.cc:431
List::getFirst
T getFirst() const
Definition: ftmpl_list.cc:279
List::length
int length() const
Definition: ftmpl_list.cc:273
List::append
void append(const T &)
Definition: ftmpl_list.cc:256
List::insert
void insert(const T &)
Definition: ftmpl_list.cc:193
Variable
factory's class for variables
Definition: factory.h:127
bigintmat
Matrices of numbers.
Definition: bigintmat.h:51
intvec
Definition: intvec.h:23
intvec::length
int length() const
Definition: intvec.h:94
ip_smatrix
Definition: matpol.h:15
coeffs.h
Coefficient rings, fields and other domains suitable for Singular polynomials.
n_Copy
static FORCE_INLINE number n_Copy(number n, const coeffs r)
return a copy of 'n'
Definition: coeffs.h:448
n_convFactoryNSingN
static FORCE_INLINE number n_convFactoryNSingN(const CanonicalForm n, const coeffs r)
Definition: coeffs.h:975
n_Invers
static FORCE_INLINE number n_Invers(number a, const coeffs r)
return the multiplicative inverse of 'a'; raise an error if 'a' is not invertible
Definition: coeffs.h:561
n_GreaterZero
static FORCE_INLINE BOOLEAN n_GreaterZero(number n, const coeffs r)
ordered fields: TRUE iff 'n' is positive; in Z/pZ: TRUE iff 0 < m <= roundedBelow(p/2),...
Definition: coeffs.h:491
n_Div
static FORCE_INLINE number n_Div(number a, number b, const coeffs r)
return the quotient of 'a' and 'b', i.e., a/b; raises an error if 'b' is not invertible in r exceptio...
Definition: coeffs.h:612
n_convSingNFactoryN
static FORCE_INLINE CanonicalForm n_convSingNFactoryN(number n, BOOLEAN setChar, const coeffs r)
Definition: coeffs.h:978
n_IsZero
static FORCE_INLINE BOOLEAN n_IsZero(number n, const coeffs r)
TRUE iff 'n' represents the zero element.
Definition: coeffs.h:461
n_Size
static FORCE_INLINE int n_Size(number n, const coeffs r)
return a non-negative measure for the complexity of n; return 0 only when n represents zero; (used fo...
Definition: coeffs.h:567
n_Delete
static FORCE_INLINE void n_Delete(number *p, const coeffs r)
delete 'p'
Definition: coeffs.h:452
n_Init
static FORCE_INLINE number n_Init(long i, const coeffs r)
a number representing i in the given coeff field/ring r
Definition: coeffs.h:535
n_IsOne
static FORCE_INLINE BOOLEAN n_IsOne(number n, const coeffs r)
TRUE iff 'n' represents the one element.
Definition: coeffs.h:465
nCoeff_is_transExt
static FORCE_INLINE BOOLEAN nCoeff_is_transExt(const coeffs r)
TRUE iff r represents a transcendental extension field.
Definition: coeffs.h:915
pa
static BOOLEAN pa(leftv res, leftv args)
Definition: cohomo.cc:3723
pb
static BOOLEAN pb(leftv res, leftv args)
Definition: cohomo.cc:3747
Warn
#define Warn
Definition: emacs.cc:77
WarnS
#define WarnS
Definition: emacs.cc:78
iter
CFFListIterator iter
Definition: facAbsBiFact.cc:53
alpha
Variable alpha
Definition: facAbsBiFact.cc:51
absFactorize
CFAFList absFactorize(const CanonicalForm &G)
absolute factorization of a multivariate poly over Q
Definition: facAbsFact.cc:262
s
const CanonicalForm int s
Definition: facAbsFact.cc:51
res
CanonicalForm res
Definition: facAbsFact.cc:60
w
const CanonicalForm & w
Definition: facAbsFact.cc:51
mipo
CanonicalForm mipo
Definition: facAlgExt.cc:57
v
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:39
j
int j
Definition: facHensel.cc:110
rootOf
Variable FACTORY_PUBLIC rootOf(const CanonicalForm &, char name='@')
returns a symbolic root of polynomial with name name Use it to define algebraic variables
Definition: variable.cc:162
prune
void FACTORY_PUBLIC prune(Variable &alpha)
Definition: variable.cc:261
errorreported
VAR short errorreported
Definition: feFopen.cc:23
WerrorS
void WerrorS(const char *s)
Definition: feFopen.cc:24
nfMinPoly
STATIC_VAR int nfMinPoly[16]
Definition: ffields.cc:545
flint_mpoly.h
flintconv.h
This file is work in progress and currently not part of the official Singular.
VAR
#define VAR
Definition: globaldefs.h:5
idIs0
BOOLEAN idIs0(ideal h)
returns true if h is the zero ideal
Definition: simpleideals.cc:917
ivCopy
intvec * ivCopy(const intvec *o)
Definition: intvec.h:145
IMATELEM
#define IMATELEM(M, I, J)
Definition: intvec.h:85
G
STATIC_VAR TreeM * G
Definition: janet.cc:31
h
STATIC_VAR Poly * h
Definition: janet.cc:971
Q
STATIC_VAR jList * Q
Definition: janet.cc:30
mpNew
matrix mpNew(int r, int c)
create a r x c zero-matrix
Definition: matpol.cc:37
iiWriteMatrix
void iiWriteMatrix(matrix im, const char *n, int dim, const ring r, int spaces)
set spaces to zero by default
Definition: matpol.cc:834
MATELEM
#define MATELEM(mat, i, j)
1-based access to matrix
Definition: matpol.h:29
assume
#define assume(x)
Definition: mod2.h:389
p_GetComp
#define p_GetComp(p, r)
Definition: monomials.h:64
pIter
#define pIter(p)
Definition: monomials.h:37
pNext
#define pNext(p)
Definition: monomials.h:36
p_GetCoeff
#define p_GetCoeff(p, r)
Definition: monomials.h:50
pGetCoeff
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy
Definition: monomials.h:44
LinTree::init
void init()
Definition: lintree.cc:864
coeffs
The main handler for Singular numbers which are suitable for Singular polynomials.
ndConvSingNFactoryN
CanonicalForm ndConvSingNFactoryN(number, BOOLEAN, const coeffs)
Definition: numbers.cc:313
omAlloc0
#define omAlloc0(size)
Definition: omAllocDecl.h:211
NULL
#define NULL
Definition: omList.c:12
p_IsPurePower
int p_IsPurePower(const poly p, const ring r)
return i, if head depends only on var(i)
Definition: p_polys.cc:1226
p_Cleardenom_n
void p_Cleardenom_n(poly ph, const ring r, number &c)
Definition: p_polys.cc:2950
p_GcdMon
poly p_GcdMon(poly f, poly g, const ring r)
polynomial gcd for f=mon
Definition: p_polys.cc:4879
p_Content
void p_Content(poly ph, const ring r)
Definition: p_polys.cc:2291
p_Norm
void p_Norm(poly p1, const ring r)
Definition: p_polys.cc:3715
p_Div_mm
poly p_Div_mm(poly p, const poly m, const ring r)
divide polynomial by monomial
Definition: p_polys.cc:1534
p_Cleardenom
poly p_Cleardenom(poly p, const ring r)
Definition: p_polys.cc:2841
p_One
poly p_One(const ring r)
Definition: p_polys.cc:1313
p_Sub
poly p_Sub(poly p1, poly p2, const ring r)
Definition: p_polys.cc:1986
p_NSet
poly p_NSet(number n, const ring r)
returns the poly representing the number n, destroys n
Definition: p_polys.cc:1469
p_polys.h
p_Neg
static poly p_Neg(poly p, const ring r)
Definition: p_polys.h:1105
pLength
static int pLength(poly a)
Definition: p_polys.h:188
p_Add_q
static poly p_Add_q(poly p, poly q, const ring r)
Definition: p_polys.h:934
p_Mult_q
static poly p_Mult_q(poly p, poly q, const ring r)
Definition: p_polys.h:1112
p_Write
void p_Write(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:342
p_SetExp
static unsigned long p_SetExp(poly p, const unsigned long e, const unsigned long iBitmask, const int VarOffset)
set a single variable exponent @Note: VarOffset encodes the position in p->exp
Definition: p_polys.h:486
p_SetComp
static unsigned long p_SetComp(poly p, unsigned long c, ring r)
Definition: p_polys.h:245
p_Setm
static void p_Setm(poly p, const ring r)
Definition: p_polys.h:231
p_SetCoeff
static number p_SetCoeff(poly p, number n, ring r)
Definition: p_polys.h:410
p_GetExp
static long p_GetExp(const poly p, const unsigned long iBitmask, const int VarOffset)
get a single variable exponent @Note: the integer VarOffset encodes:
Definition: p_polys.h:467
p_IsConstant
static BOOLEAN p_IsConstant(const poly p, const ring r)
Definition: p_polys.h:1962
p_Delete
static void p_Delete(poly *p, const ring r)
Definition: p_polys.h:899
p_Write0
void p_Write0(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:332
pp_Mult_qq
static poly pp_Mult_qq(poly p, poly q, const ring r)
Definition: p_polys.h:1149
p_Copy
static poly p_Copy(poly p, const ring r)
returns a copy of p
Definition: p_polys.h:844
p_Totaldegree
static long p_Totaldegree(poly p, const ring r)
Definition: p_polys.h:1505
p_Test
#define p_Test(p, r)
Definition: p_polys.h:159
__p_Mult_nn
#define __p_Mult_nn(p, n, r)
Definition: p_polys.h:969
p_wrp
void p_wrp(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:373
currRing
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:13
pTotaldegree
static long pTotaldegree(poly p)
Definition: polys.h:282
pWrite
void pWrite(poly p)
Definition: polys.h:308
mark
std::pair< int, int > mark
Definition: ppinitialReduction.h:7
NUM
@ NUM
Definition: readcf.cc:170
StringSetS
void StringSetS(const char *st)
Definition: reporter.cc:128
StringAppendS
void StringAppendS(const char *st)
Definition: reporter.cc:107
feNotImplemented
const char feNotImplemented[]
Definition: reporter.cc:54
PrintS
void PrintS(const char *s)
Definition: reporter.cc:284
StringEndS
char * StringEndS()
Definition: reporter.cc:151
PrintLn
void PrintLn()
Definition: reporter.cc:310
Werror
void Werror(const char *fmt,...)
Definition: reporter.cc:189
rChar
int rChar(ring r)
Definition: ring.cc:713
rField_is_Zp_a
static BOOLEAN rField_is_Zp_a(const ring r)
Definition: ring.h:529
rField_is_Z
static BOOLEAN rField_is_Z(const ring r)
Definition: ring.h:509
rField_is_Zp
static BOOLEAN rField_is_Zp(const ring r)
Definition: ring.h:500
rField_is_Zn
static BOOLEAN rField_is_Zn(const ring r)
Definition: ring.h:512
rPar
static int rPar(const ring r)
(r->cf->P)
Definition: ring.h:599
rInternalChar
static int rInternalChar(const ring r)
Definition: ring.h:689
rField_is_Q_a
static BOOLEAN rField_is_Q_a(const ring r)
Definition: ring.h:539
rField_is_Q
static BOOLEAN rField_is_Q(const ring r)
Definition: ring.h:506
rField_is_GF
static BOOLEAN rField_is_GF(const ring r)
Definition: ring.h:521
rVar
static short rVar(const ring r)
#define rVar(r) (r->N)
Definition: ring.h:592
count
int status int void size_t count
Definition: si_signals.h:59
idInit
ideal idInit(int idsize, int rank)
initialise an ideal / module
Definition: simpleideals.cc:35
id_Delete
void id_Delete(ideal *h, ring r)
deletes an ideal/module/matrix
Definition: simpleideals.cc:123
idShow
void idShow(const ideal id, const ring lmRing, const ring tailRing, const int debugPrint)
Definition: simpleideals.cc:57
idSkipZeroes
void idSkipZeroes(ideal ide)
gives an ideal/module the minimal possible size
Definition: simpleideals.cc:181
IDELEMS
#define IDELEMS(i)
Definition: simpleideals.h:23
id_Test
#define id_Test(A, lR)
Definition: simpleideals.h:87
R
#define R
Definition: sirandom.c:27
M
#define M
Definition: sirandom.c:25
loop
#define loop
Definition: structs.h:75
nf
Definition: gnumpfl.cc:25
gcd
int gcd(int a, int b)
Definition: walkSupport.cc:836
Generated on Thu Apr 27 2023 08:16:02 for My Project by   doxygen 1.9.4 for Singular