[5a66d0] | 1 | #include <stdio.h> |
---|
| 2 | #include <stdlib.h> |
---|
| 3 | #include <string.h> |
---|
| 4 | #include <time.h> |
---|
| 5 | |
---|
[762407] | 6 | #include "config.h" |
---|
[b1dfaf] | 7 | #include <kernel/mod2.h> |
---|
| 8 | #include <omalloc/omalloc.h> |
---|
[5ff1d3] | 9 | |
---|
[0fb34ba] | 10 | #include <coeffs/numbers.h> |
---|
[5ff1d3] | 11 | |
---|
[737a68] | 12 | #include <kernel/polys.h> |
---|
[0fb34ba] | 13 | #include <polys/monomials/ring.h> |
---|
[fe66ba8] | 14 | #include <polys/monomials/p_polys.h> |
---|
[5ff1d3] | 15 | #include <polys/kbuckets.h> |
---|
| 16 | |
---|
| 17 | #include <kernel/ideals.h> |
---|
| 18 | #include <kernel/longrat.h> |
---|
| 19 | #include <kernel/febase.h> |
---|
| 20 | #include <kernel/kutil.h> |
---|
| 21 | |
---|
| 22 | #include "subexpr.h" |
---|
| 23 | |
---|
| 24 | |
---|
| 25 | #include "janet.h" |
---|
[5a66d0] | 26 | |
---|
| 27 | #if (defined(__CYGWIN__)) |
---|
| 28 | #include <ctype.h> |
---|
| 29 | #endif |
---|
[5ff1d3] | 30 | |
---|
[5a66d0] | 31 | #include <stdarg.h> |
---|
| 32 | |
---|
| 33 | |
---|
| 34 | //------GLOBALS------- |
---|
| 35 | static int m_s,v_s,vectorized,VarN1,offset; |
---|
| 36 | static jList *T,*Q; |
---|
| 37 | static TreeM *G; |
---|
| 38 | static Poly *phD; |
---|
| 39 | static NodeM *FreeNodes; |
---|
| 40 | static int degree_compatible; |
---|
| 41 | static int (*ListGreatMove)(jList *,jList *,poly); |
---|
| 42 | static int Mask[8]={0x80,0x40,0x20,0x10,0x8,0x4,0x2,0x1}; |
---|
| 43 | |
---|
| 44 | //#define DebugPrint |
---|
| 45 | |
---|
| 46 | //#define pow_(x) pTotaldegree((x)) |
---|
[31f1850] | 47 | //#define pow_(x) p_Deg((x,currRing)) |
---|
[5a66d0] | 48 | pFDegProc jDeg; |
---|
[08db0a] | 49 | #define pow_(x) jDeg((x),currRing) |
---|
[5a66d0] | 50 | |
---|
| 51 | void Debug() |
---|
| 52 | { |
---|
| 53 | LCI it=T->root; |
---|
| 54 | |
---|
[ca4d0e6] | 55 | PrintS("T==================================\n"); |
---|
[5a66d0] | 56 | while (it) |
---|
| 57 | { |
---|
| 58 | pWrite(it->info->root); |
---|
| 59 | it=it->next; |
---|
| 60 | } |
---|
| 61 | |
---|
| 62 | it=Q->root; |
---|
| 63 | |
---|
[ca4d0e6] | 64 | PrintS("Q==================================\n"); |
---|
[5a66d0] | 65 | while (it) |
---|
| 66 | { |
---|
| 67 | if (it->info->root) pWrite(it->info->root); |
---|
| 68 | else |
---|
| 69 | { |
---|
| 70 | Print("%d.........",it->info->prolonged); |
---|
| 71 | pWrite(it->info->history); |
---|
| 72 | } |
---|
| 73 | it=it->next; |
---|
| 74 | } |
---|
[ca4d0e6] | 75 | PrintS("===================================\n"); |
---|
[5a66d0] | 76 | } |
---|
| 77 | |
---|
| 78 | int ReducePolyLead(Poly *x,Poly *y) |
---|
| 79 | { |
---|
| 80 | if (!x->root || !y->root) |
---|
| 81 | return 0; |
---|
| 82 | |
---|
| 83 | /* poly b1=pDivide(x->root,y->root); |
---|
| 84 | |
---|
| 85 | number gcd=nGcd(pGetCoeff(x->root),pGetCoeff(y->root),currRing); |
---|
| 86 | |
---|
| 87 | number a1=nDiv(pGetCoeff(y->root),gcd); |
---|
| 88 | pGetCoeff(b1)=nDiv(pGetCoeff(x->root),gcd); |
---|
| 89 | |
---|
| 90 | x->root=pMult_nn(x->root,a1); |
---|
| 91 | nDelete(&a1); |
---|
| 92 | |
---|
| 93 | x->root=pMinus_mm_Mult_qq(x->root,b1,y->root); |
---|
| 94 | |
---|
| 95 | pDelete(&b1); |
---|
| 96 | */ |
---|
| 97 | #if 1 |
---|
| 98 | if (x->root_b==NULL) |
---|
| 99 | { |
---|
| 100 | if (x->root_l<=0) x->root_l=pLength(x->root); |
---|
| 101 | x->root_b=kBucketCreate(currRing); |
---|
| 102 | kBucketInit(x->root_b,x->root,x->root_l); |
---|
| 103 | } |
---|
| 104 | number coef; |
---|
| 105 | if (y->root_l<=0) y->root_l=pLength(y->root); |
---|
| 106 | coef=kBucketPolyRed(x->root_b,y->root,y->root_l,NULL); |
---|
| 107 | nDelete(&coef); |
---|
| 108 | x->root=kBucketGetLm(x->root_b); |
---|
| 109 | if (x->root==NULL) |
---|
| 110 | { |
---|
| 111 | kBucketDestroy(&x->root_b); |
---|
| 112 | x->root_b=NULL; |
---|
| 113 | x->root_l=0; |
---|
| 114 | } |
---|
| 115 | #else |
---|
| 116 | x->root=ksOldSpolyRed(y->root,x->root,NULL); |
---|
| 117 | #endif |
---|
[6e0750] | 118 | // if (x->root) p_Content(x->root,currRing); |
---|
[5a66d0] | 119 | // if (x->root) pSimpleContent(x->root,5); |
---|
| 120 | |
---|
| 121 | return 1; |
---|
| 122 | } |
---|
| 123 | |
---|
| 124 | int ReducePoly(Poly *x,poly from,Poly *y) |
---|
| 125 | { |
---|
| 126 | if (!x->root || !y->root) |
---|
| 127 | return 0; |
---|
| 128 | |
---|
| 129 | /* poly b1=pDivide(from,y->root); |
---|
| 130 | |
---|
| 131 | number gcd=nGcd(pGetCoeff(from),pGetCoeff(y->root),currRing); |
---|
| 132 | |
---|
| 133 | number a1=nDiv(pGetCoeff(y->root),gcd); |
---|
| 134 | pGetCoeff(b1)=nDiv(pGetCoeff(from),gcd); |
---|
| 135 | |
---|
| 136 | x->root=pMult_nn(x->root,a1); |
---|
| 137 | nDelete(&a1);*/ |
---|
| 138 | |
---|
| 139 | // x->root=pMinus_mm_Mult_qq(x->root,b1,y->root); |
---|
| 140 | // pDelete(&b1); |
---|
| 141 | |
---|
| 142 | ksOldSpolyTail(y->root,x->root,from,NULL,currRing); |
---|
| 143 | y->root_l=0; |
---|
| 144 | |
---|
| 145 | return 1; |
---|
| 146 | } |
---|
| 147 | |
---|
| 148 | void PNF(Poly *p, TreeM *F) |
---|
| 149 | { |
---|
| 150 | if (p->root==NULL) return; |
---|
| 151 | |
---|
| 152 | Poly *f; |
---|
| 153 | BOOLEAN done=FALSE; |
---|
| 154 | poly temp=p->root; |
---|
| 155 | |
---|
| 156 | // if (TEST_OPT_PROT) { PrintS("r"); mflush(); } |
---|
| 157 | int count=0; |
---|
| 158 | poly pp=p->root; |
---|
| 159 | int old_size=nSize(pGetCoeff(pp)); |
---|
| 160 | p->root_l=0; |
---|
| 161 | while(temp->next) |
---|
| 162 | { |
---|
| 163 | f=is_div_(F,temp->next); |
---|
| 164 | if (f) |
---|
| 165 | { |
---|
| 166 | if (ReducePoly(p,temp,f)) //temp->next |
---|
| 167 | { |
---|
| 168 | count++; |
---|
| 169 | //if (TEST_OPT_PROT) { PrintS("-"); mflush(); } |
---|
| 170 | if ((f!=NULL) |
---|
| 171 | && (count>20) |
---|
| 172 | && (nSize(pGetCoeff(pp))>old_size) |
---|
| 173 | ) |
---|
| 174 | { |
---|
| 175 | //pSimpleContent(pp,2); |
---|
[6e0750] | 176 | p_Content(pp,currRing); |
---|
[5a66d0] | 177 | count=0; |
---|
| 178 | // old_size=nSize(pGetCoeff(pp)); |
---|
| 179 | } |
---|
| 180 | } |
---|
| 181 | done=TRUE; |
---|
| 182 | } |
---|
| 183 | else |
---|
| 184 | temp=temp->next; |
---|
| 185 | } |
---|
| 186 | |
---|
[6e0750] | 187 | if (done) p_Content(p->root,currRing); |
---|
[5a66d0] | 188 | //if (done) pSimpleContent(p->root,-1); |
---|
| 189 | pTest(p->root); |
---|
| 190 | } |
---|
| 191 | |
---|
| 192 | void NFL(Poly *p, TreeM *F) |
---|
| 193 | { |
---|
| 194 | Poly *f; |
---|
| 195 | int g1,f1,gg; |
---|
| 196 | |
---|
| 197 | if ((f=is_div_(F,p->lead))==NULL) return; |
---|
| 198 | |
---|
| 199 | int pX=pow_(p->lead); |
---|
| 200 | int phX=pow_(p->history); |
---|
| 201 | |
---|
| 202 | if (pX!=phX) |
---|
| 203 | { |
---|
| 204 | int phF=pow_(f->history); |
---|
| 205 | if (pX >= (phX+phF)) |
---|
| 206 | { |
---|
| 207 | pDelete(&p->root); |
---|
| 208 | //p->root=NULL; |
---|
| 209 | return; |
---|
| 210 | } |
---|
| 211 | |
---|
| 212 | /* poly p2=pInit(); |
---|
| 213 | pLcm(p->history,f->history,p2); |
---|
| 214 | pSetm(p2); |
---|
| 215 | |
---|
| 216 | if (pLmCmp(p->root,p2) > 0) |
---|
| 217 | { |
---|
| 218 | pLmDelete(&p2); |
---|
| 219 | pDelete(&p->root); |
---|
| 220 | //p->root=NULL; |
---|
| 221 | return; |
---|
| 222 | } |
---|
| 223 | |
---|
| 224 | pLmDelete(&p2); |
---|
| 225 | */ |
---|
| 226 | /* for(int i=0, gg=0 ; i<currRing->N;i++) |
---|
| 227 | if ((g1=pGetExp(p->history,i+1)) > (f1=pGetExp(f->history,i+1))) |
---|
| 228 | gg+=g1; |
---|
| 229 | else gg+=f1; |
---|
| 230 | |
---|
| 231 | if (pX > gg) |
---|
| 232 | { |
---|
| 233 | pDelete(&p->root); |
---|
| 234 | //x->root=NULL; |
---|
| 235 | return; |
---|
| 236 | } |
---|
| 237 | */ |
---|
| 238 | int pF=pow_(f->lead); |
---|
| 239 | |
---|
| 240 | if ((pX == pF) && (pF == phF)) |
---|
| 241 | { |
---|
| 242 | pLmDelete(&f->history); |
---|
| 243 | f->history=pCopy(p->history); |
---|
| 244 | } |
---|
| 245 | } |
---|
| 246 | |
---|
| 247 | //if (TEST_OPT_PROT) { PrintS("R"); mflush(); } |
---|
| 248 | int old_size, count; |
---|
| 249 | count=0; |
---|
| 250 | while(f && p->root) |
---|
| 251 | { |
---|
[ca4d0e6] | 252 | // PrintS("R"); |
---|
[5a66d0] | 253 | // if (TEST_OPT_PROT) { PrintS("R"); mflush(); } |
---|
| 254 | #if 0 |
---|
| 255 | old_size=nSize(pGetCoeff(p->root)); |
---|
| 256 | #endif |
---|
| 257 | if (ReducePolyLead(p,f) == 0) break; |
---|
| 258 | if (p->root!=NULL) |
---|
| 259 | { |
---|
| 260 | count++; |
---|
| 261 | #if 0 |
---|
| 262 | if ((count>4) && (3<nSize(pGetCoeff(p->root))) |
---|
| 263 | && (nSize(pGetCoeff(p->root))>old_size)) |
---|
| 264 | { |
---|
| 265 | pSimpleContent(p->root,old_size); |
---|
| 266 | count=0; |
---|
| 267 | } |
---|
| 268 | #else |
---|
| 269 | if (count>500) |
---|
| 270 | { |
---|
| 271 | kBucketClear(p->root_b,&p->root,&p->root_l); |
---|
[fe66ba8] | 272 | p_SimpleContent(p->root,2,currRing); |
---|
[5a66d0] | 273 | kBucketInit(p->root_b,p->root,p->root_l); |
---|
| 274 | count=0; |
---|
[ca4d0e6] | 275 | //PrintS("."); |
---|
[5a66d0] | 276 | } |
---|
| 277 | #endif |
---|
| 278 | f=is_div_(F,p->root); |
---|
| 279 | } |
---|
| 280 | } |
---|
| 281 | #if 1 |
---|
| 282 | if (p->root_b!=NULL) |
---|
| 283 | { |
---|
| 284 | kBucketClear(p->root_b,&p->root,&p->root_l); |
---|
| 285 | kBucketDestroy(&p->root_b); |
---|
| 286 | p->root_b=NULL; |
---|
| 287 | } |
---|
| 288 | #endif |
---|
| 289 | |
---|
| 290 | if (!p->root) |
---|
| 291 | return; |
---|
| 292 | |
---|
| 293 | InitHistory(p); |
---|
| 294 | InitProl(p); |
---|
| 295 | InitLead(p); |
---|
| 296 | p->changed=1; |
---|
| 297 | |
---|
[6e0750] | 298 | p_Content(p->root,currRing); |
---|
[5a66d0] | 299 | //pSimpleContent(p->root,-1); |
---|
| 300 | pTest(p->root); |
---|
| 301 | } |
---|
| 302 | |
---|
| 303 | int ValidatePoly(Poly *x, TreeM *F) |
---|
| 304 | { |
---|
| 305 | Poly *f,*g; |
---|
| 306 | int g1,f1; |
---|
| 307 | |
---|
| 308 | if (x->root) return 1; |
---|
| 309 | |
---|
| 310 | g=is_present(T,x->history); //it's a prolongation - do we have a parent ? |
---|
| 311 | |
---|
| 312 | if (!g) return 0; //if not - kill him ! |
---|
| 313 | |
---|
| 314 | poly lmX=pDivide(x->lead,g->root); |
---|
| 315 | pGetCoeff(lmX)=nInit(1); |
---|
| 316 | |
---|
| 317 | /* if ((f=is_div_(F,lmX)) != NULL) |
---|
| 318 | { |
---|
| 319 | int pX=pow_(lmX); |
---|
| 320 | int phX=pow_(x->history); |
---|
| 321 | |
---|
| 322 | if (pX!=phX) |
---|
| 323 | { |
---|
| 324 | int phF=pow_(f->history); |
---|
| 325 | if (pX >= (phX+phF)) |
---|
| 326 | { |
---|
| 327 | pLmDelete(&lmX); |
---|
| 328 | //x->root=NULL; |
---|
| 329 | return 0; |
---|
| 330 | } |
---|
| 331 | |
---|
| 332 | for(int i=0, gg=0 ; i<currRing->N;i++) |
---|
| 333 | if ((g1=pGetExp(x->history,i+1)) > (f1=pGetExp(f->history,i+1))) |
---|
| 334 | gg+=g1; |
---|
| 335 | else |
---|
| 336 | gg+=f1; |
---|
| 337 | |
---|
| 338 | if (pX > gg) |
---|
| 339 | { |
---|
| 340 | pLmDelete(&lmX); |
---|
| 341 | return 0; |
---|
| 342 | } |
---|
| 343 | int pF=pow_(f->root); |
---|
| 344 | |
---|
| 345 | if ((pX == pF) && (pF == phF)) |
---|
| 346 | f->history=x->history; |
---|
| 347 | } |
---|
| 348 | } |
---|
| 349 | |
---|
| 350 | pLmDelete(&lmX); |
---|
| 351 | |
---|
| 352 | */ |
---|
| 353 | x->root=pCopy(g->root); |
---|
| 354 | x->root_l=g->root_l; |
---|
| 355 | |
---|
| 356 | x->root=pMult(x->root,lmX); |
---|
| 357 | |
---|
| 358 | pTest(x->root); |
---|
| 359 | |
---|
| 360 | x->prolonged=-1; |
---|
| 361 | |
---|
| 362 | return 1; |
---|
| 363 | } |
---|
| 364 | |
---|
| 365 | Poly *NewPoly(poly p) |
---|
| 366 | { |
---|
| 367 | Poly *beg=(Poly *)GCM(sizeof(Poly)); |
---|
| 368 | |
---|
| 369 | beg->root=p;//(p == NULL ? pInit() : p); |
---|
| 370 | beg->root_b=NULL; |
---|
| 371 | beg->root_l=0; |
---|
| 372 | beg->history=NULL;//pInit(); |
---|
| 373 | beg->lead=NULL; |
---|
| 374 | beg->mult=(char *)GCMA(sizeof(char)*2*offset); |
---|
| 375 | |
---|
| 376 | for (int i=0; i < currRing->N; i++) |
---|
| 377 | { |
---|
| 378 | ClearMult(beg,i); |
---|
| 379 | ClearProl(beg,i); |
---|
| 380 | }; |
---|
| 381 | |
---|
| 382 | beg->prolonged=-1; |
---|
| 383 | |
---|
| 384 | return beg; |
---|
| 385 | } |
---|
| 386 | |
---|
| 387 | void DestroyPoly(Poly *x) |
---|
| 388 | { |
---|
| 389 | pDelete(&x->root); |
---|
| 390 | pDelete(&x->history); |
---|
| 391 | if (x->lead) pDelete(&x->lead); |
---|
| 392 | GCF(x->mult); |
---|
| 393 | GCF(x); |
---|
| 394 | } |
---|
| 395 | |
---|
| 396 | void ControlProlong(Poly *x) |
---|
| 397 | { |
---|
| 398 | for (int i = 0; i< offset; i++) |
---|
| 399 | { |
---|
| 400 | (x->mult+offset)[i]&=~((x->mult)[i]); |
---|
| 401 | // if (!GetMult(x,i) && !GetProl(x,i)) |
---|
| 402 | // ProlVar(x,i); |
---|
| 403 | } |
---|
| 404 | } |
---|
| 405 | |
---|
| 406 | void InitHistory(Poly *p) |
---|
| 407 | { |
---|
| 408 | if (p->history) pLmDelete(&p->history); |
---|
| 409 | p->history=pLmInit(p->root); |
---|
| 410 | p->changed=0; |
---|
| 411 | } |
---|
| 412 | |
---|
| 413 | void InitLead(Poly *p) |
---|
| 414 | { |
---|
| 415 | if (p->lead) pLmDelete(&p->lead); |
---|
| 416 | p->lead=pLmInit(p->root); |
---|
| 417 | p->prolonged=-1; |
---|
| 418 | } |
---|
| 419 | |
---|
| 420 | void InitProl(Poly *p) |
---|
| 421 | { |
---|
| 422 | memset(p->mult+offset,0,sizeof(char)*offset); |
---|
| 423 | } |
---|
| 424 | |
---|
| 425 | int GetMult(Poly *x,int i) |
---|
| 426 | { |
---|
| 427 | return x->mult[i/8] & Mask[i%8]; |
---|
| 428 | } |
---|
| 429 | |
---|
| 430 | void SetMult(Poly *x,int i) |
---|
| 431 | { |
---|
| 432 | x->mult[i/8] |= Mask[i%8]; |
---|
| 433 | } |
---|
| 434 | |
---|
| 435 | void ClearMult(Poly *x,int i) |
---|
| 436 | { |
---|
| 437 | x->mult[i/8] &= ~Mask[i%8]; |
---|
| 438 | } |
---|
| 439 | |
---|
| 440 | int GetProl(Poly *x, int i) |
---|
| 441 | { |
---|
| 442 | return (x->mult+offset)[i/8] & Mask[i%8]; |
---|
| 443 | } |
---|
| 444 | |
---|
| 445 | void SetProl(Poly *x, int i) |
---|
| 446 | { |
---|
| 447 | (x->mult+offset)[i/8] |= Mask[i%8]; |
---|
| 448 | } |
---|
| 449 | |
---|
| 450 | void ClearProl(Poly *x, int i) |
---|
| 451 | { |
---|
| 452 | (x->mult+offset)[i/8] &= ~Mask[i%8]; |
---|
| 453 | } |
---|
| 454 | |
---|
| 455 | int LengthCompare(poly p1,poly p2) |
---|
| 456 | { |
---|
| 457 | do |
---|
| 458 | { |
---|
| 459 | if (p1 == NULL) return 1; |
---|
| 460 | if (p2 == NULL) return 0; |
---|
| 461 | pIter(p1); |
---|
| 462 | pIter(p2); |
---|
| 463 | }while(p1 && p2); |
---|
| 464 | return 1; |
---|
| 465 | } |
---|
| 466 | |
---|
| 467 | int ProlCompare(Poly *item1, Poly *item2) |
---|
| 468 | { |
---|
| 469 | switch(pLmCmp(item1->lead,item2->lead)) |
---|
| 470 | { |
---|
| 471 | case -1: |
---|
| 472 | return 1; |
---|
| 473 | |
---|
| 474 | case 1: |
---|
| 475 | return 0; |
---|
| 476 | |
---|
| 477 | default: |
---|
| 478 | if ((item1->root_l<=0)||(item2->root_l<=0)) |
---|
| 479 | return LengthCompare(item1->root,item2->root); |
---|
| 480 | return item1->root_l<=item2->root_l; |
---|
| 481 | } |
---|
| 482 | } |
---|
| 483 | |
---|
| 484 | void ProlVar(Poly *temp,int i) |
---|
| 485 | { |
---|
| 486 | Poly *Pr; |
---|
| 487 | |
---|
| 488 | if (!GetProl(temp,i) && !GetMult(temp,i)) |
---|
| 489 | { |
---|
| 490 | Pr=NewPoly(); |
---|
| 491 | SetProl(temp,i); |
---|
| 492 | |
---|
| 493 | Pr->prolonged=i; |
---|
| 494 | Pr->history=pLmInit(temp->history); |
---|
| 495 | Pr->lead=pLmInit(temp->lead); |
---|
| 496 | pIncrExp(Pr->lead,i+1); |
---|
| 497 | pSetm(Pr->lead); |
---|
| 498 | InitProl(temp); |
---|
| 499 | |
---|
| 500 | Pr->changed=0; |
---|
| 501 | // pTest(Pr->root); |
---|
| 502 | InsertInCount(Q,Pr); |
---|
| 503 | } |
---|
| 504 | } |
---|
| 505 | |
---|
| 506 | void DestroyListNode(ListNode *x) |
---|
| 507 | { |
---|
| 508 | DestroyPoly(x->info); |
---|
| 509 | GCF(x); |
---|
| 510 | } |
---|
| 511 | |
---|
| 512 | ListNode* CreateListNode(Poly *x) |
---|
| 513 | { |
---|
| 514 | ListNode* ret=(ListNode *)GCM(sizeof(ListNode)); |
---|
| 515 | ret->info=x; |
---|
| 516 | ret->next=NULL; |
---|
| 517 | return ret; |
---|
| 518 | } |
---|
| 519 | |
---|
| 520 | |
---|
| 521 | Poly *FindMinList(jList *L) |
---|
| 522 | { |
---|
| 523 | LI min=&(L->root); |
---|
| 524 | LI l; |
---|
| 525 | LCI xl; |
---|
| 526 | Poly *x; |
---|
| 527 | |
---|
| 528 | if (degree_compatible) |
---|
| 529 | { |
---|
| 530 | while ((*min) && ((*min)->info->root == NULL)) |
---|
| 531 | min=&((*min)->next); |
---|
| 532 | } |
---|
| 533 | |
---|
| 534 | if (!(*min)) return NULL; |
---|
| 535 | |
---|
| 536 | l=&((*min)->next); |
---|
| 537 | |
---|
| 538 | while (*l) |
---|
| 539 | { |
---|
| 540 | if ((*l)->info->root != NULL) |
---|
| 541 | { |
---|
| 542 | if (ProlCompare((*l)->info,(*min)->info)) |
---|
| 543 | min=l; |
---|
| 544 | } |
---|
| 545 | |
---|
| 546 | l=&((*l)->next); |
---|
| 547 | } |
---|
| 548 | x=(*min)->info; |
---|
| 549 | xl=*min; |
---|
| 550 | *min=(*min)->next; |
---|
| 551 | GCF(xl); |
---|
| 552 | |
---|
| 553 | return x; |
---|
| 554 | } |
---|
| 555 | |
---|
| 556 | void InsertInList(jList *x,Poly *y) |
---|
| 557 | { |
---|
| 558 | ListNode *ins; |
---|
| 559 | LI ix=&(x->root); |
---|
| 560 | |
---|
| 561 | while (*ix) |
---|
| 562 | { |
---|
| 563 | if (pLmCmp(y->lead,(*ix)->info->lead) == -1) |
---|
| 564 | ix=(ListNode **)&((*ix)->next); |
---|
| 565 | else |
---|
| 566 | break; |
---|
| 567 | } |
---|
| 568 | |
---|
| 569 | ins=CreateListNode(y); |
---|
| 570 | ins->next=(ListNode *)(*ix); |
---|
| 571 | *ix=ins; |
---|
| 572 | return; |
---|
| 573 | } |
---|
| 574 | |
---|
| 575 | void InsertInCount(jList *x,Poly *y) |
---|
| 576 | { |
---|
| 577 | ListNode *ins; |
---|
| 578 | LI ix=&(x->root); |
---|
| 579 | |
---|
| 580 | ins=CreateListNode(y); |
---|
| 581 | ins->next=(ListNode *)(*ix); |
---|
| 582 | *ix=ins; |
---|
| 583 | return; |
---|
| 584 | } |
---|
| 585 | |
---|
| 586 | int ListGreatMoveOrder(jList *A,jList *B,poly x) |
---|
| 587 | { |
---|
| 588 | LCI y=A->root; |
---|
| 589 | |
---|
| 590 | if (!y || pLmCmp(y->info->lead,x) < 0) return 0; |
---|
| 591 | |
---|
| 592 | while(y && pLmCmp(y->info->lead,x) >= 0) |
---|
| 593 | { |
---|
| 594 | InsertInCount(B,y->info); |
---|
| 595 | A->root=y->next; |
---|
| 596 | GCF(y); |
---|
| 597 | y=A->root; |
---|
| 598 | } |
---|
| 599 | |
---|
| 600 | return 1; |
---|
| 601 | } |
---|
| 602 | |
---|
| 603 | int ListGreatMoveDegree(jList *A,jList *B,poly x) |
---|
| 604 | { |
---|
| 605 | LCI y=A->root; |
---|
| 606 | int pow_x=pow_(x); |
---|
| 607 | |
---|
| 608 | if (!y || pow_(y->info->lead) <= pow_x) return 0; |
---|
| 609 | |
---|
| 610 | while(y && pow_(y->info->lead) > pow_x) |
---|
| 611 | { |
---|
| 612 | InsertInCount(B,y->info); |
---|
| 613 | A->root=y->next; |
---|
| 614 | GCF(y); |
---|
| 615 | y=A->root; |
---|
| 616 | } |
---|
| 617 | |
---|
| 618 | return 1; |
---|
| 619 | } |
---|
| 620 | |
---|
| 621 | int CountList(jList *Q) |
---|
| 622 | { |
---|
| 623 | int i=0; |
---|
| 624 | LCI y=Q->root; |
---|
| 625 | |
---|
| 626 | while(y) |
---|
| 627 | { |
---|
| 628 | i++; |
---|
| 629 | y=y->next; |
---|
| 630 | } |
---|
| 631 | |
---|
| 632 | return i; |
---|
| 633 | } |
---|
| 634 | |
---|
| 635 | void NFListQ() |
---|
| 636 | { |
---|
| 637 | LCI ll; |
---|
| 638 | int p,p1; |
---|
| 639 | LI l; |
---|
| 640 | |
---|
| 641 | do |
---|
| 642 | { |
---|
| 643 | if (!Q->root) break; |
---|
| 644 | |
---|
| 645 | ll=Q->root; |
---|
| 646 | |
---|
| 647 | p=pow_(Q->root->info->lead); |
---|
| 648 | |
---|
| 649 | while (ll) |
---|
| 650 | { |
---|
| 651 | int ploc=pow_(ll->info->lead); |
---|
| 652 | if (ploc < p) p=ploc; |
---|
| 653 | ll=ll->next; |
---|
| 654 | } |
---|
| 655 | |
---|
| 656 | p1=1; |
---|
| 657 | |
---|
| 658 | l=&(Q->root); |
---|
| 659 | |
---|
| 660 | while (*l) |
---|
| 661 | { |
---|
[ca4d0e6] | 662 | // PrintS("*"); |
---|
[5a66d0] | 663 | int ploc=pow_((*l)->info->lead); |
---|
| 664 | |
---|
| 665 | if (ploc == p) |
---|
| 666 | { |
---|
| 667 | if (!ValidatePoly((*l)->info,G)) |
---|
| 668 | { |
---|
| 669 | ll=(*l); |
---|
| 670 | *l=(*l)->next; |
---|
| 671 | DestroyListNode(ll); |
---|
| 672 | continue; |
---|
| 673 | }; |
---|
| 674 | |
---|
| 675 | (*l)->info->changed=0; |
---|
[ca4d0e6] | 676 | // PrintS("!"); |
---|
[5a66d0] | 677 | NFL((*l)->info,G); |
---|
[ca4d0e6] | 678 | // PrintS("$"); |
---|
[5a66d0] | 679 | if (!(*l)->info->root) |
---|
| 680 | { |
---|
| 681 | ll=(*l); |
---|
| 682 | *l=(*l)->next; |
---|
| 683 | DestroyListNode(ll); |
---|
| 684 | continue; |
---|
| 685 | }; |
---|
| 686 | p1=0; |
---|
| 687 | } |
---|
| 688 | |
---|
| 689 | l=&((*l)->next); |
---|
| 690 | } |
---|
| 691 | }while(p1); |
---|
[ca4d0e6] | 692 | // PrintLn(); |
---|
[5a66d0] | 693 | } |
---|
| 694 | |
---|
| 695 | |
---|
| 696 | void ForEachPNF(jList *x,int i) |
---|
| 697 | { |
---|
| 698 | LCI y=x->root; |
---|
| 699 | |
---|
| 700 | while(y) |
---|
| 701 | { |
---|
| 702 | if (pow_(y->info->root) == i) PNF(y->info,G); |
---|
| 703 | y=y->next; |
---|
| 704 | } |
---|
| 705 | } |
---|
| 706 | |
---|
| 707 | void ForEachControlProlong(jList *x) |
---|
| 708 | { |
---|
| 709 | LCI y=x->root; |
---|
| 710 | |
---|
| 711 | while(y) |
---|
| 712 | { |
---|
| 713 | ControlProlong(y->info); |
---|
| 714 | y=y->next; |
---|
| 715 | } |
---|
| 716 | } |
---|
| 717 | |
---|
| 718 | void DestroyList(jList *x) |
---|
| 719 | { |
---|
| 720 | LCI y=x->root,z; |
---|
| 721 | |
---|
| 722 | while(y) |
---|
| 723 | { |
---|
| 724 | z=y->next; |
---|
| 725 | DestroyPoly(y->info); |
---|
| 726 | GCF(y); |
---|
| 727 | y=z; |
---|
| 728 | } |
---|
| 729 | |
---|
| 730 | GCF(x); |
---|
| 731 | } |
---|
| 732 | |
---|
| 733 | Poly* is_present(jList *F,poly x) |
---|
| 734 | { |
---|
| 735 | LCI iF=F->root; |
---|
| 736 | while(iF) |
---|
| 737 | if (pLmCmp(iF->info->root,x) == 0) |
---|
| 738 | return iF->info; |
---|
| 739 | else iF=iF->next; |
---|
| 740 | |
---|
| 741 | return NULL; |
---|
| 742 | } |
---|
| 743 | |
---|
| 744 | int GB_length() |
---|
| 745 | { |
---|
| 746 | LCI iT=T->root; |
---|
| 747 | int l=0; |
---|
| 748 | |
---|
| 749 | while(iT) |
---|
| 750 | { |
---|
| 751 | if (pow_(iT->info->lead) == pow_(iT->info->history)) |
---|
| 752 | ++l; |
---|
| 753 | iT=iT->next; |
---|
| 754 | } |
---|
| 755 | |
---|
| 756 | return l; |
---|
| 757 | } |
---|
| 758 | |
---|
| 759 | static Poly *temp_l; |
---|
| 760 | |
---|
| 761 | NodeM* create() |
---|
| 762 | { |
---|
| 763 | NodeM *y; |
---|
| 764 | |
---|
| 765 | if (FreeNodes == NULL) |
---|
| 766 | { |
---|
| 767 | y=(NodeM *)GCM(sizeof(NodeM)); |
---|
| 768 | } |
---|
| 769 | else |
---|
| 770 | { |
---|
| 771 | y=FreeNodes; |
---|
| 772 | FreeNodes=FreeNodes->left; |
---|
| 773 | } |
---|
| 774 | |
---|
| 775 | y->left=y->right=NULL; |
---|
| 776 | y->ended=NULL; |
---|
| 777 | return y; |
---|
| 778 | } |
---|
| 779 | |
---|
| 780 | void DestroyFreeNodes() |
---|
| 781 | { |
---|
| 782 | NodeM *y; |
---|
| 783 | |
---|
| 784 | while((y=FreeNodes)!=NULL) |
---|
| 785 | { |
---|
| 786 | FreeNodes=FreeNodes->left; |
---|
| 787 | GCF(y); |
---|
| 788 | } |
---|
| 789 | } |
---|
| 790 | |
---|
| 791 | static void go_right(NodeM *current,poly_function disp) |
---|
| 792 | { |
---|
| 793 | if (current) |
---|
| 794 | { |
---|
| 795 | go_right(current->left,disp); |
---|
| 796 | if (current->ended) disp(current->ended); |
---|
| 797 | go_right(current->right,disp); |
---|
| 798 | } |
---|
| 799 | } |
---|
| 800 | |
---|
| 801 | void ForEach(TreeM *t,poly_function disp) |
---|
| 802 | { |
---|
| 803 | go_right(t->root,disp); |
---|
| 804 | } |
---|
| 805 | |
---|
| 806 | void DestroyTree(NodeM *G) |
---|
| 807 | { |
---|
| 808 | if (G) |
---|
| 809 | { |
---|
| 810 | DestroyTree(G->left); |
---|
| 811 | DestroyTree(G->right); |
---|
| 812 | G->left=FreeNodes; |
---|
| 813 | FreeNodes=G; |
---|
| 814 | } |
---|
| 815 | } |
---|
| 816 | |
---|
| 817 | void Define(TreeM **G) |
---|
| 818 | { |
---|
| 819 | *G=(TreeM *)GCM(sizeof(TreeM)); |
---|
| 820 | (*G)->root=create(); |
---|
| 821 | } |
---|
| 822 | |
---|
| 823 | int sp_div(poly m1,poly m2,int from) |
---|
| 824 | { |
---|
| 825 | |
---|
| 826 | if (pow_(m2) == 0 && pow_(m1)) return 0; |
---|
| 827 | |
---|
| 828 | for(int k=from; k < currRing->N; k++) |
---|
| 829 | if (pGetExp(m1,k+1) < pGetExp(m2,k+1)) return 0; |
---|
| 830 | |
---|
| 831 | return 1; |
---|
| 832 | } |
---|
| 833 | |
---|
| 834 | void div_l(poly item, NodeM *x,int from) |
---|
| 835 | { |
---|
| 836 | if (x && !temp_l) |
---|
| 837 | { |
---|
| 838 | div_l(item,x->left,from); |
---|
| 839 | if ((x->ended) && sp_div(item,x->ended->root,from)) |
---|
| 840 | { |
---|
| 841 | temp_l=x->ended; |
---|
| 842 | return; |
---|
| 843 | }; |
---|
| 844 | div_l(item,x->right,from); |
---|
| 845 | } |
---|
| 846 | } |
---|
| 847 | |
---|
| 848 | Poly* is_div_upper(poly item, NodeM *x,int from) |
---|
| 849 | { |
---|
| 850 | temp_l=NULL; |
---|
| 851 | div_l(item,x,from); |
---|
| 852 | return temp_l; |
---|
| 853 | } |
---|
| 854 | |
---|
| 855 | Poly* is_div_(TreeM *tree, poly item) |
---|
| 856 | { |
---|
| 857 | int power_tmp,i,i_con=currRing->N-1; |
---|
| 858 | NodeM *curr=tree->root; |
---|
| 859 | |
---|
| 860 | if (!curr) return NULL; |
---|
| 861 | if (pow_(item) == 0) return NULL; |
---|
| 862 | |
---|
| 863 | for ( ; i_con>=0 && !pGetExp(item,i_con+1) ; i_con--) |
---|
| 864 | ; |
---|
| 865 | |
---|
| 866 | for (i=0; i <= i_con ; i++) |
---|
| 867 | { |
---|
| 868 | power_tmp=pGetExp(item,i+1); |
---|
| 869 | |
---|
| 870 | while (power_tmp) |
---|
| 871 | { |
---|
| 872 | if (curr->ended) return curr->ended; |
---|
| 873 | |
---|
| 874 | if (!curr->left) |
---|
| 875 | { |
---|
| 876 | if (curr->right) |
---|
| 877 | return is_div_upper(item,curr->right,i); //?????? |
---|
| 878 | return NULL; |
---|
| 879 | } |
---|
| 880 | |
---|
| 881 | curr=curr->left; |
---|
| 882 | power_tmp--; |
---|
| 883 | } |
---|
| 884 | |
---|
| 885 | if (curr->ended) return curr->ended; |
---|
| 886 | |
---|
| 887 | if (!curr->right) return NULL; |
---|
| 888 | |
---|
| 889 | curr=curr->right; |
---|
| 890 | } |
---|
| 891 | |
---|
| 892 | if (curr->ended) return curr->ended; |
---|
| 893 | else return NULL; |
---|
| 894 | } |
---|
| 895 | |
---|
| 896 | static void ClearMultiplicative(NodeM *xx,int i) |
---|
| 897 | { |
---|
| 898 | if (!xx) return; |
---|
| 899 | |
---|
| 900 | while (xx->left) |
---|
| 901 | { |
---|
| 902 | ClearMultiplicative(xx->right, i); |
---|
| 903 | xx = xx->left; |
---|
| 904 | } |
---|
| 905 | if ((xx->ended) && (GetMult(xx->ended,i))) |
---|
| 906 | { |
---|
| 907 | ClearMult(xx->ended,i); |
---|
| 908 | ProlVar(xx->ended,i); |
---|
| 909 | } |
---|
| 910 | else |
---|
| 911 | ClearMultiplicative(xx->right,i); |
---|
| 912 | } |
---|
| 913 | //====================================================== |
---|
| 914 | void insert_(TreeM **tree, Poly *item) |
---|
| 915 | { |
---|
| 916 | int power_tmp,i,i_con=currRing->N-1; |
---|
| 917 | NodeM *curr=(*tree)->root; |
---|
| 918 | |
---|
| 919 | for ( ; (i_con>=0) && !pGetExp(item->root,i_con+1) ; i_con--) |
---|
| 920 | SetMult(item,i_con); |
---|
| 921 | |
---|
| 922 | for (i = 0; i<= i_con; i++) |
---|
| 923 | //<= |
---|
| 924 | { |
---|
| 925 | power_tmp=pGetExp(item->root,i+1); |
---|
| 926 | |
---|
| 927 | ClearMult(item,i); |
---|
| 928 | |
---|
| 929 | while (power_tmp) |
---|
| 930 | { |
---|
| 931 | if (!curr->left) |
---|
| 932 | { |
---|
| 933 | SetMult(item,i); |
---|
| 934 | ClearMultiplicative(curr->right,i); |
---|
| 935 | curr->left=create(); |
---|
| 936 | }; |
---|
| 937 | curr=curr->left; |
---|
| 938 | power_tmp--; |
---|
| 939 | }; |
---|
| 940 | |
---|
| 941 | if (i<i_con) |
---|
| 942 | { |
---|
| 943 | if (!curr->left) SetMult(item,i); |
---|
| 944 | if (!curr->right) curr->right=create(); |
---|
| 945 | curr=curr->right; |
---|
| 946 | |
---|
| 947 | ProlVar(item,i); |
---|
| 948 | } |
---|
| 949 | } |
---|
| 950 | |
---|
| 951 | curr->ended=item; |
---|
| 952 | } |
---|
| 953 | |
---|
| 954 | void Initialization(char *Ord) |
---|
| 955 | { |
---|
| 956 | offset=(currRing->N % 8 == 0) ? (currRing->N/8)*8 : (currRing->N/8+1)*8; |
---|
| 957 | if (strstr(Ord,"dp\0") || strstr(Ord,"Dp\0")) |
---|
| 958 | { |
---|
| 959 | degree_compatible=1; |
---|
[31f1850] | 960 | jDeg=p_Deg; |
---|
[5a66d0] | 961 | ListGreatMove=ListGreatMoveDegree; |
---|
| 962 | } |
---|
| 963 | else |
---|
| 964 | { |
---|
| 965 | degree_compatible=0; |
---|
[99bdcf] | 966 | jDeg=p_Totaldegree; |
---|
[5a66d0] | 967 | ListGreatMove=ListGreatMoveOrder; |
---|
| 968 | } |
---|
| 969 | |
---|
| 970 | Define(&G); |
---|
[e8a9f3] | 971 | } |
---|
[5a66d0] | 972 | |
---|
| 973 | static Poly *h,*f; |
---|
| 974 | |
---|
| 975 | void insert_in_G(Poly *x) |
---|
| 976 | { |
---|
| 977 | insert_(&G,x); |
---|
| 978 | } |
---|
| 979 | |
---|
| 980 | void T2G(); |
---|
| 981 | |
---|
| 982 | void Q2TG() |
---|
| 983 | { |
---|
| 984 | LCI t; |
---|
| 985 | Poly *x; |
---|
| 986 | |
---|
| 987 | while (Q->root) |
---|
| 988 | { |
---|
| 989 | t=Q->root; |
---|
| 990 | x=t->info; |
---|
| 991 | insert_(&G,x); |
---|
| 992 | InsertInList(T,x); |
---|
| 993 | Q->root=t->next; |
---|
| 994 | GCF(t); |
---|
| 995 | } |
---|
| 996 | } |
---|
| 997 | |
---|
[eb72ba1] | 998 | int ComputeBasis(jList *_lT,jList *_lQ) |
---|
[5a66d0] | 999 | { |
---|
| 1000 | int gb_l,i,ret_value=1; |
---|
| 1001 | |
---|
[eb72ba1] | 1002 | T=_lT; Q=_lQ; |
---|
[5a66d0] | 1003 | |
---|
| 1004 | // Debug(); |
---|
| 1005 | |
---|
| 1006 | while((h=FindMinList(Q))!=NULL) |
---|
| 1007 | { |
---|
[ca4d0e6] | 1008 | // PrintS("New element\n"); |
---|
[5a66d0] | 1009 | // Debug(); |
---|
| 1010 | |
---|
| 1011 | if (!degree_compatible) |
---|
| 1012 | { |
---|
| 1013 | if (!ValidatePoly(h,G)) |
---|
| 1014 | { |
---|
| 1015 | DestroyPoly(h); |
---|
| 1016 | continue; |
---|
| 1017 | } |
---|
| 1018 | |
---|
| 1019 | h->changed=0; |
---|
| 1020 | |
---|
| 1021 | NFL(h,G); |
---|
| 1022 | |
---|
| 1023 | if (!h->root) |
---|
| 1024 | { |
---|
| 1025 | DestroyPoly(h); |
---|
| 1026 | continue; |
---|
| 1027 | } |
---|
| 1028 | } |
---|
| 1029 | |
---|
| 1030 | if (h->root) |
---|
| 1031 | { |
---|
| 1032 | if (pIsConstant(h->root)) |
---|
| 1033 | { |
---|
| 1034 | WarnS("Constant in basis\n"); |
---|
| 1035 | return 0; |
---|
| 1036 | } |
---|
| 1037 | |
---|
| 1038 | if (h->changed && ListGreatMove(T,Q,h->root)) |
---|
| 1039 | { |
---|
[ca4d0e6] | 1040 | // PrintS("<-\n"); |
---|
[5a66d0] | 1041 | DestroyTree(G->root); |
---|
| 1042 | G->root=create(); |
---|
| 1043 | T2G(); |
---|
| 1044 | } |
---|
| 1045 | } |
---|
| 1046 | |
---|
[ca4d0e6] | 1047 | // PrintS("PNF\n"); |
---|
[5a66d0] | 1048 | PNF(h,G); |
---|
| 1049 | // Print("{%d}\n",pow_(h->root)); |
---|
| 1050 | insert_(&G,h); |
---|
| 1051 | InsertInList(T,h); |
---|
| 1052 | |
---|
[ca4d0e6] | 1053 | // PrintS("For each PNF\n"); |
---|
[5a66d0] | 1054 | if (degree_compatible) |
---|
| 1055 | ForEachPNF(T,pow_(h->root)); |
---|
| 1056 | |
---|
[ca4d0e6] | 1057 | // PrintS("Control of prolongations\n"); |
---|
[5a66d0] | 1058 | if (h->changed) |
---|
| 1059 | ForEachControlProlong(T); |
---|
| 1060 | else |
---|
| 1061 | ControlProlong(h); |
---|
| 1062 | |
---|
| 1063 | // Debug(); |
---|
| 1064 | |
---|
[ca4d0e6] | 1065 | // PrintS("NFListQ\n"); |
---|
[5a66d0] | 1066 | if (degree_compatible) |
---|
| 1067 | NFListQ(); |
---|
| 1068 | //Debug(); |
---|
| 1069 | } |
---|
| 1070 | |
---|
| 1071 | // gb_l=GB_length(); |
---|
| 1072 | |
---|
| 1073 | Print("Length of Janet basis: %d\n",CountList(T)); |
---|
| 1074 | // Print("Length of Groebner basis: %d\n",gb_l); |
---|
| 1075 | |
---|
| 1076 | DestroyTree(G->root); |
---|
| 1077 | GCF(G); |
---|
| 1078 | DestroyFreeNodes(); |
---|
| 1079 | |
---|
| 1080 | return 1; |
---|
| 1081 | } |
---|
| 1082 | |
---|
| 1083 | void T2G() |
---|
| 1084 | { |
---|
| 1085 | LCI i=T->root; |
---|
| 1086 | while (i) |
---|
| 1087 | { |
---|
| 1088 | insert_(&G,i->info); |
---|
| 1089 | i=i->next; |
---|
| 1090 | } |
---|
| 1091 | } |
---|