My Project
Loading...
Searching...
No Matches
cohomo.cc
Go to the documentation of this file.
1/****************************************
2* Computer Algebra System SINGULAR *
3****************************************/
4/*
5* ABSTRACT - Stainly
6*/
7
8#include "kernel/mod2.h"
9
10#if !defined(__CYGWIN__) || defined(STATIC_VERSION)
11// acces from a module to routines from the main program
12// does not work on windows (restrict of the dynamic linker),
13// a static version is required:
14// ./configure --with-builtinmodules=cohomo,...
15
16
17#include "omalloc/omalloc.h"
18#include "misc/mylimits.h"
20#include <assert.h>
21#include <unistd.h>
22
26#include "cohomo.h"//for my thing
27#include "kernel/GBEngine/tgb.h"//
28#include "Singular/ipid.h"//for ggetid
31#include "polys/simpleideals.h"
32#include "Singular/lists.h"
33#include "kernel/linear_algebra/linearAlgebra.h"//for printNumber
34#include "kernel/GBEngine/kstd1.h"//for gb
35#include <kernel/ideals.h>
36#if 1
37
41#include <coeffs/numbers.h>
42#include <vector>
43#include <Singular/ipshell.h>
45#include <time.h>
46
47/***************************print(only for debugging)***********************************************/
48//print vector of integers.
49static void listprint(std::vector<int> vec)
50{
51 unsigned i;
52 for(i=0;i<vec.size();i++)
53 {
54 Print(" _[%d]=%d\n",i+1,vec[i]);
55 PrintLn();
56 }
57 if(vec.size()==0)
58 {
59 PrintS(" _[1]= \n");
60 PrintLn();
61 }
62}
63
64//print vector of vectors of integers.
65static void listsprint(std::vector<std::vector<int> > posMat)
66{
67 unsigned i;
68 for(i=0;i<posMat.size();i++)
69 {
70 Print("[%d]:\n",i+1);
71 listprint(posMat[i]);
72 Print("\n");
73 PrintLn();
74 }
75 if(posMat.size()==0)
76 {
77 PrintS("[1]:\n");
78 PrintLn();
79 }
80}
81
82
83//print ideal.
84static void id_print(ideal h)
85{
86 int i;
87 for(i=0;i<IDELEMS(h);i++)
88 {
89 Print(" [%d]\n",i+1);
90 pWrite(h->m[i]);
91 PrintLn();
92 }
93}
94
95//only for T^2,
96//print vector of polynomials.
97static void lpprint( std::vector<poly> pv)
98{
99 for(unsigned i=0;i<pv.size();i++)
100 {
101 Print(" _[%d]=",i+1);
102 pWrite(pv[i]);
103 }
104 if(pv.size()==0)
105 {
106 PrintS(" _[1]= \n");
107 PrintLn();
108 }
109}
110
111//print vector of vectors of polynomials.
112static void lpsprint(std::vector<std::vector<poly> > pvs)
113{
114 for(unsigned i=0;i<pvs.size();i++)
115 {
116 Print("[%d]:\n",i+1);
117 lpprint(pvs[i]);
118 Print("\n");
119 PrintLn();
120 }
121 if(pvs.size()==0)
122 {
123 PrintS("[1]:\n");
124 PrintLn();
125 }
126}
127
128/*************operations for vectors (regard vectors as sets)*********/
129
130//returns true if integer n is in vector vec,
131//otherwise, returns false
132static bool IsinL(int a, std::vector<int> vec)
133{
134 unsigned i;
135 for(i=0;i<vec.size();i++)
136 {
137 if(a==vec[i])
138 {
139 return true;
140 }
141 }
142 return false;
143}
144
145//returns intersection of vectors p and q,
146//returns empty if they are disjoint
147static std::vector<int> vecIntersection(std::vector<int> p, std::vector<int> q)
148{
149 unsigned i;
150 std::vector<int> inte;
151 for(i=0;i<p.size();i++)
152 {
153 if(IsinL(p[i],q))
154 inte.push_back(p[i]);
155 }
156 return inte;
157}
158
159//returns true if vec1 is contained in vec2
160static bool vsubset(std::vector<int> vec1, std::vector<int> vec2)
161{
162 int i;
163 if(vec1.size()>vec2.size())
164 return false;
165 for(i=0;i<vec1.size();i++)
166 {
167 if(!IsinL(vec1[i],vec2))
168 return false;
169 }
170 return true;
171}
172
173//not strictly equal(order doesn't matter)
174static bool vEvl(std::vector<int> vec1, std::vector<int> vec2)
175{
176 if(vec1.size()==0 && vec2.size()==0)
177 return true;
178 if(vsubset(vec1,vec2)&&vsubset(vec2,vec1))
179 return true;
180 return false;
181}
182
183//the length of vec must be same to it of the elements of vecs
184//returns true if vec is as same as some element of vecs ((not strictly same))
185//returns false if vec is not in vecs
186static bool vInvsl(std::vector<int> vec, std::vector<std::vector<int> > vecs)
187{
188 int i;
189 for(i=0;i<vecs.size();i++)
190 {
191 if(vEvl(vec,vecs[i]))
192 {
193 return true;
194 }
195 }
196 return false;
197}
198
199//returns the union of two vectors(as the union of sets)
200static std::vector<int> vecUnion(std::vector<int> vec1, std::vector<int> vec2)
201{
202 std::vector<int> vec=vec1;
203 unsigned i;
204 for(i=0;i<vec2.size();i++)
205 {
206 if(!IsinL(vec2[i],vec))
207 vec.push_back(vec2[i]);
208 }
209 return vec;
210}
211
212static std::vector<int> vecMinus(std::vector<int> vec1,std::vector<int> vec2)
213{
214 std::vector<int> vec;
215 for(unsigned i=0;i<vec1.size();i++)
216 {
217 if(!IsinL(vec1[i],vec2))
218 {
219 vec.push_back(vec1[i]);
220 }
221 }
222 return vec;
223}
224
225static std::vector<std::vector<int> > vsMinusv(std::vector<std::vector<int> > vecs, std::vector<int> vec)
226{
227 int i;
228 std::vector<std::vector<int> > rem;
229 for(i=0;i<vecs.size();i++)
230 {
231 if(!vEvl(vecs[i],vec))
232 {
233 rem.push_back(vecs[i]);
234 }
235 }
236 return (rem);
237}
238
239static std::vector<std::vector<int> > vsUnion(std::vector<std::vector<int> > vs1, std::vector<std::vector<int> > vs2)
240{
241 int i;
242 std::vector<std::vector<int> > vs=vs1;
243 for(i=0;i<vs2.size();i++)
244 {
245 if(!vInvsl(vs2[i],vs))
246 {
247 vs.push_back(vs2[i]);
248 }
249 }
250 return vs;
251}
252
253static std::vector<std::vector<int> > vsIntersection(std::vector<std::vector<int> > vs1, std::vector<std::vector<int> > vs2)
254{
255 int i;
256 std::vector<std::vector<int> > vs;
257 for(i=0;i<vs2.size();i++)
258 {
259 if(vInvsl(vs2[i],vs1))
260 {
261 vs.push_back(vs2[i]);
262 }
263 }
264 return vs;
265}
266
267/*************************************for transition between ideal and vectors******************************************/
268
269//P should be monomial,
270// vector version of poly support(poly p)
271static std::vector<int> support1(poly p)
272{
273 int j;
274 std::vector<int> supset;
275 if(p==0) return supset;
276 for(j=1;j<=rVar(currRing);j++)
277 {
278 if(pGetExp(p,j)>0)
279 {
280 supset.push_back(j);
281 }
282 }
283 return (supset);
284}
285
286//simplicial complex(the faces set is ideal h)
287static std::vector<std::vector<int> > supports(ideal h)
288{
289 std::vector<std::vector<int> > vecs;
290 std::vector<int> vec;
291 if(!idIs0(h))
292 {
293 for(int s=0;s<IDELEMS(h);s++)
294 {
295 vec=support1(h->m[s]);
296 vecs.push_back(vec);
297 }
298 }
299 return vecs;
300}
301
302// only for eqsolve1
303// p could be any polynomial
304static std::vector<int> support2(poly p)
305{
306 int j;
307 poly q;
308 std::vector<int> supset;
309 for(j=1;j<=rVar(currRing);j++)
310 {
311 q=pCopy(p);
312 while (q!=NULL)
313 {
314 if(p_GetExp(q,j,currRing)!=0)
315 {
316 supset.push_back(j);
317 break;
318 }
319 q=pNext(q);
320 }
321 }
322 return (supset);
323}
324
325//the supports of ideal
326static std::vector<std::vector<int> > supports2(ideal h)
327{
328 std::vector<std::vector<int> > vecs;
329 std::vector<int> vec;
330 if(!idIs0(h))
331 {
332 for(int s=0;s<IDELEMS(h);s++)
333 {
334 vec=support2(h->m[s]);
335 vecs.push_back(vec);
336 }
337 }
338 return vecs;
339}
340
341//convert the vector(vbase[i] are the coefficients of x_{i+1}) to a polynomial w.r.t. current ring
342//vector vbase has length of currRing->N.
343static poly pMake(std::vector<int> vbase)
344{
345 int n=vbase.size(); poly p,q=0;
346 for(int i=0;i<n;i++)
347 {
348 if(vbase[i]!=0)
349 {
350 p = pOne();pSetExp(p, i+1, 1);pSetm(p);pSetCoeff(p, nInit(vbase[i]));
351 q = pAdd(q, p);
352 }
353 }
354 return q;
355}
356
357//convert the vectors to a ideal(for T^1)
358static ideal idMake(std::vector<std::vector<int> > vecs)
359{
360 int lv=vecs.size(), i;
361 poly p;
362 ideal id_re=idInit(1,1);
363 for(i=0;i<lv;i++)
364 {
365 p=pMake(vecs[i]);
366 idInsertPoly(id_re, p);
367 }
368 idSkipZeroes(id_re);
369 return id_re;
370}
371
372/*****************************quotient ring of two ideals*********************/
373
374//the quotient ring of h1 respect to h2
375static ideal idmodulo(ideal h1,ideal h2)
376{
377 ideal gb=kStd(h2,NULL,testHomog,NULL,NULL,0,0,NULL);
378 idSkipZeroes(gb);
379 ideal idq=kNF(gb,NULL,h1);
380 idSkipZeroes(idq);
381 return idq;
382}
383
384//returns the coeff of the monomial of polynomial p which involves the mth varialbe
385//assume the polynomial p has form of y1+y2+...
386static int pcoef(poly p, int m)
387{
388 int i,co; poly q=pCopy(p);
389 for(i=1;i<=currRing->N;i++)
390 {
391 if(p_GetExp(q,m,currRing)!=0)
392 {
393 co=n_Int(pGetCoeff(q),currRing->cf);
394 return co;
395 }
396 else
397 q=pNext(q);
398 }
399 if(q!=NULL)
400 co=0;
401 return co;
402}
403
404//returns true if p involves the mth variable of the current ring
405static bool vInp(int m,poly p)
406{
407 poly q=pCopy(p);
408 while (q!=NULL)
409 {
410 if(p_GetExp(q,m,currRing)!=0)
411 {
412 return true;
413 }
414 q=pNext(q);
415 }
416 return false;
417}
418
419//returns the vector w.r.t. polynomial p
420static std::vector<int> vMake(poly p)
421{
422 int i;
423 std::vector<int> vbase;
424 for(i=1;i<=currRing->N;i++)
425 {
426 if(vInp(i,p))
427 {
428 vbase.push_back(pcoef(p,i));
429 }
430 else
431 {
432 vbase.push_back(0);
433 }
434 }
435 return (vbase);
436}
437
438//returns the vectors w.r.t. ideal h
439static std::vector<std::vector<int> > vsMake(ideal h)
440{
441 std::vector<int> vec;
442 std::vector<std::vector<int> > vecs;
443 int i;
444 for(i=0;i<IDELEMS(h);i++)
445 {
446 vec=vMake(h->m[i]);
447 vecs.push_back(vec);
448 }
449 return vecs;
450}
451
452//the quotient ring of two ideals which are represented by vectors,
453//the result is also represented by vector.
454static std::vector<std::vector<int> > vecqring(std::vector<std::vector<int> > vec1, std::vector<std::vector<int> > vec2)
455{
456 ideal h1=idMake(vec1), h2=idMake(vec2);
457 ideal h=idmodulo(h1,h2);
458 std::vector<std::vector<int> > vecs= vsMake(h);
459 return vecs;
460}
461
462/****************************************************************/
463//construct a monomial based on the support of it
464//returns a squarefree monomial
465static poly pMaken(std::vector<int> vbase)
466{
467 int n=vbase.size();
468 poly p,q=pOne();
469 for(int i=0;i<n;i++)
470 {
471 p = pOne();pSetExp(p, vbase[i], 1);pSetm(p);pSetCoeff(p, nInit(1));
472 //pWrite(p);
473 q=pp_Mult_mm(q,p,currRing);
474 }
475 return q;
476}
477
478// returns a ideal according to a set of supports
479static ideal idMaken(std::vector<std::vector<int> > vecs)
480{
481 ideal id_re=idInit(1,1);
482 poly p;
483 int i,lv=vecs.size();
484 for(i=0;i<lv;i++)
485 {
486 p=pMaken(vecs[i]);
487 idInsertPoly(id_re, p);
488 }
489 idSkipZeroes(id_re);
490 //id_print(id_re);
491 return id_re;
492}
493
494/********************************new version for stanley reisner ideal ***********************************************/
495
496static std::vector<std::vector<int> > b_subsets(std::vector<int> vec)
497{
498 int i,j;
499 std::vector<int> bv;
500 std::vector<std::vector<int> > vecs;
501 for(i=0;i<vec.size();i++)
502 {
503 bv.push_back(vec[i]);
504 vecs.push_back(bv);
505 bv.clear();
506 }
507 //listsprint(vecs);
508 for(i=0;i<vecs.size();i++)
509 {
510 for(j=i+1;j<vecs.size();j++)
511 {
512 bv=vecUnion(vecs[i], vecs[j]);
513 if(!vInvsl(bv,vecs))
514 vecs.push_back(bv);
515 }
516 }
517 //listsprint(vecs);
518 return(vecs);
519}
520
521//the number of the variables
522static int idvert(ideal h)
523{
524 int i, j, vert=0;
525 if(idIs0(h))
526 return vert;
527 for(i=currRing->N;i>0;i--)
528 {
529 for(j=0;j<IDELEMS(h);j++)
530 {
531 if(pGetExp(h->m[j],i)>0)
532 {
533 vert=i;
534 return vert;
535 }
536 }
537 }
538 return vert;
539}
540
541static int pvert(poly p)
542{
543 int i, vert=0;
544 for(i=currRing->N;i>0;i--)
545 {
546 if(pGetExp(p,i)>0)
547 {
548 vert=i;
549 return vert;
550 }
551 }
552 return vert;
553}
554
555/*
556//full complex
557static std::vector<std::vector<int> > fullcomplex(ideal h)
558{
559 int vert=vertnum(h), i, j;
560 //Print("there are %d vertices\n", vert);
561 std::vector<std::vector<int> > fmons;
562 std::vector<int> pre;
563 for(i=1;i<=vert;i++)
564 {
565 pre.push_back(i);
566 }
567 fmons=b_subsets(pre);
568 return fmons;
569}*/
570
571/*
572//all the squarefree monomials whose degree is less or equal to n
573static std::vector<std::vector<int> > sfrmons(ideal h, int n)
574{
575 int vert=vertnum(h), i, j, time=0;
576 std::vector<std::vector<int> > fmons, pres, pres0, pres1;
577 std::vector<int> pre;
578 for(i=1;i<=vert;i++)
579 {
580 pre.push_back(i);
581 pres0.push_back(pre);
582 }
583 pres=pres0;
584 for(i=0;i<size(pres),time<=n;i++)
585 {
586 time++;
587 pre=pres[i];
588 for(j=0;j<size(pres0);j++)
589 {
590 pre=vecUnion(pre, pres0[j]);
591 if(pre.)
592 }
593 }
594 return fmons;
595}
596*/
597
598/*
599static ideal id_complement(ideal h)
600{
601 int i,j;
602 std::vector<std::vector<int> > full=fullcomplex(h), hvs=supports(h), res;
603 for(i=0;i<full.size();i++)
604 {
605 if(!vInvsl(full[i], hvs))
606 {
607 res.push_back(full[i]);
608 }
609 }
610 return idMaken(res);
611}*/
612
613
614/*****************About simplicial complex and stanley-reisner ideal and ring **************************/
615
616//h1 minus h2
617static ideal idMinus(ideal h1,ideal h2)
618{
619 ideal h=idInit(1,1);
620 int i,j,eq=0;
621 for(i=0;i<IDELEMS(h1);i++)
622 {
623 eq=0;
624 for(j=0;j<IDELEMS(h2);j++)
625 {
626 if(p_EqualPolys(pCopy(h1->m[i]),pCopy(h2->m[j]), currRing))
627 {
628 eq=1;
629 break;
630 }
631 }
632 if(eq==0)
633 {
634 idInsertPoly(h, pCopy(h1->m[i]));
635 }
636 }
638 return h;
639}
640
641//If poly P is squarefree, returns 1
642//returns 0 otherwise,
643static bool p_Ifsfree(poly P)
644{
645 int i,sf=1;
646 for(i=1;i<=rVar(currRing);i++)
647 {
648 if (pGetExp(P,i)>1)
649 {
650 sf=0;
651 break;
652 }
653 }
654 return sf;
655}
656
657//returns the set of all squarefree monomials of degree deg in ideal h
658static ideal sfreemon(ideal h,int deg)
659{
660 int j;
661 ideal temp;
662 temp=idInit(1,1);
663 if(!idIs0(h))
664 {
665 for(j=0;j<IDELEMS(h);j++)
666 {
667 if((p_Ifsfree(h->m[j]))&&(pTotaldegree(h->m[j])==deg))
668 {
669 idInsertPoly(temp, h->m[j]);
670 }
671 }
672 idSkipZeroes(temp);
673 }
674 return temp;
675}
676
677//full simplex represented by ideal.
678//(all the squarefree monomials over the polynomial ring)
679static ideal id_sfmon(ideal h)
680{
681 ideal asfmons,sfmons,mons;
682 int j, vert=idvert(h);
683 mons=id_MaxIdeal(1, currRing);
684 asfmons=sfreemon(mons,1);
685 for(j=2;j<=vert;j++)
686 {
687 mons=id_MaxIdeal(j, currRing);
688 sfmons=sfreemon(mons,j);
689 asfmons=id_Add(asfmons,sfmons,currRing);
690 }
691 return asfmons;
692}
693
694//if the input ideal is simplicial complex, returns the stanley-reisner ideal,
695//if the input ideal is stanley-reisner ideal, returns the monomial ideal according to simplicial complex.
696//(nonfaces and faces).
697//returns the complement of the ideal h (consisting of only squarefree polynomials)
698static ideal id_complement(ideal h)
699{
700 int j, vert=idvert(h);
701 ideal i1=id_sfmon(h);
702 ideal i3=idInit(1,1);
703 poly p;
704 for(j=0;j<IDELEMS(i1);j++)
705 {
706 p=pCopy(i1->m[j]);
707 if(pvert(p)<=vert)
708 {
709 idInsertPoly(i3, p);
710 }
711 }
712 ideal i2=idMinus(i3,h);
713 idSkipZeroes(i2);
714 return (i2);
715}
716
717//Returns true if p is one of the generators of ideal X
718//returns false otherwise
719static bool IsInX(poly p,ideal X)
720{
721 int i;
722 for(i=0;i<IDELEMS(X);i++)
723 {
724 if(pEqualPolys(p,X->m[i]))
725 {
726 //PrintS("yes\n");
727 return(true);
728 }
729 }
730 //PrintS("no\n");
731 return(false);
732}
733
734//returns the monomials in the quotient ring R/(h1+h2) which have degree deg.
735static ideal qringadd(ideal h1, ideal h2, int deg)
736{
737 ideal h,qrh;
738 h=idAdd(h1,h2);
739 qrh=scKBase(deg,h);
740 return qrh;
741}
742
743//returns the maximal degree of the monomials in ideal h
744static int id_maxdeg(ideal h)
745{
746 int i,max;
747 max=pTotaldegree(h->m[0]);
748 for(i=1;i<IDELEMS(h);i++)
749 {
750 if(pTotaldegree(h->m[i]) > max)
751 max=pTotaldegree(h->m[i]);
752 }
753 return (max);
754}
755
756//input ideal h (a squarefree monomial ideal) is the ideal associated to simplicial complex,
757//and returns the Stanley-Reisner ideal(minimal generators)
758static ideal idsrRing(ideal h)
759{
760 int i,n;
761 ideal pp,qq,rsr,ppp,hc=idCopy(h);
762 for(i=1;i<=rVar(currRing);i++)
763 {
764 pp=sfreemon(hc,i);
765 pp=scKBase(i,pp);//quotient ring (R/I_i)_i
766 if(!idIs0(pp))
767 {
768 pp=sfreemon(pp,i);
769 rsr=pp;
770 //Print("This is the first quotient generators %d:\n",i);
771 //id_print(rsr);
772 break;
773 }
774 }
775 for(n=i+1;n<=rVar(currRing);n++)
776 {
777 qq=sfreemon(hc,n);
778 pp=qringadd(qq,rsr,n);
779 ppp=sfreemon(pp,n);
780 rsr=idAdd(rsr,ppp);
781 }
782 idSkipZeroes(rsr);
783 return rsr;
784}
785
786//returns the set of all the polynomials could divide p
787static ideal SimFacset(poly p)
788{
789 int i,j,max=pTotaldegree(p);
790 ideal h1,mons,id_re=idInit(1,1);
791 for(i=1;i<max;i++)
792 {
793 mons=id_MaxIdeal(i, currRing);
794 h1=sfreemon(mons,i);
795
796 for(j=0;j<IDELEMS(h1);j++)
797 {
798 if(p_DivisibleBy(h1->m[j],p,currRing))
799 {
800 idInsertPoly(id_re, h1->m[j]);
801 }
802 }
803 }
804 idSkipZeroes(id_re);
805 return id_re;
806}
807
808static ideal idadda(ideal h1, ideal h2)
809{
810 ideal h=idInit(1,1);
811 for(int i=0;i<IDELEMS(h1);i++)
812 {
813 if(!IsInX(h1->m[i],h))
814 {
815 idInsertPoly(h, h1->m[i]);
816 }
817 }
818 for(int i=0;i<IDELEMS(h2);i++)
819 {
820 if(!IsInX(h2->m[i],h))
821 {
822 idInsertPoly(h, h2->m[i]);
823 }
824 }
826 return h;
827}
828
829//complicated version
830//(returns false if it is not a simplicial complex and print the simplex)
831//input h is need to be at least part of faces
832static ideal IsSimplex(ideal h)
833{
834 int i,max=id_maxdeg(h);
835 poly e=pOne();
836 ideal id_re, id_so=idCopy(h);
837 for(i=0;i<IDELEMS(h);i++)
838 {
839 id_re=SimFacset(h->m[i]);
840 if(!idIs0(id_re))
841 {
842 id_so=idadda(id_so, id_re);//idAdd(id_so,id_re);
843 }
844 }
845 idInsertPoly(id_so,e);
846 idSkipZeroes(id_so);
847 return (idMinus(id_so,h));
848}
849
850//input is the subset of the Stainley-Reisner ideal
851//returns the faces
852//is not used
853static ideal complementsimplex(ideal h)
854{
855 int i,j;poly p,e=pOne();
856 ideal h1=idInit(1,1), pp, h3;
857 for(i=1;i<=rVar(currRing);i++)
858 {
859 p = pOne(); pSetExp(p, i, 2); pSetm(p); pSetCoeff(p, nInit(1));
860 idInsertPoly(h1, p);
861 }
862 idSkipZeroes(h1);
863 ideal h2=idAdd(h,h1);
864 pp=scKBase(1,h2);
865 h3=idCopy(pp);
866 for(j=2;j<=rVar(currRing);j++)
867 {
868 pp=scKBase(j,h2);
869 h3=idAdd(h3,pp);
870 }
871 idInsertPoly(h3, e);
872 idSkipZeroes(h3);
873 return (h3);
874}
875
876static int dim_sim(ideal h)
877{
878 int dim=pTotaldegree(h->m[0]), i;
879 for(i=1; i<IDELEMS(h);i++)
880 {
881 if(dim<pTotaldegree(h->m[i]))
882 {
883 dim=pTotaldegree(h->m[i]);
884 }
885 }
886 return dim;
887}
888
889static int num4dim(ideal h, int n)
890{
891 int num=0;
892 for(int i=0; i<IDELEMS(h); i++)
893 {
894 if(pTotaldegree(h->m[i])==n)
895 {
896 num++;
897 }
898 }
899 return num;
900}
901
902/********************Procedures for T1(M method and N method) ***********/
903
904//h is ideal( monomial ideal) associated to simplicial complex
905//returns the all the monomials x^b (x^b must be able to divide
906//at least one monomial in Stanley-Reisner ring)
907//not so efficient
908static ideal findb(ideal h)
909{
910 ideal ib=id_sfmon(h), nonf=id_complement(h), bset=idInit(1,1);
911 poly e=pOne();
912 int i,j;
913 for(i=0;i<IDELEMS(ib);i++)
914 {
915 for(j=0;j<IDELEMS(nonf);j++)
916 {
917 if(p_DivisibleBy(ib->m[i],nonf->m[j],currRing))
918 {
919 idInsertPoly(bset, ib->m[i]);
920 break;
921 }
922 }
923 }
924 idInsertPoly(bset,e);
925 idSkipZeroes(bset);
926 return bset;
927}
928
929//h is ideal(monomial ideal associated to simplicial complex
930//1.poly S is x^b
931//2.and deg(x^a)=deg(x^b)
932//3.x^a and x^a have disjoint supports
933//returns all the possible x^a according conditions 1. 2. 3.
934static ideal finda(ideal h,poly S,int ddeg)
935{
936 poly e=pOne();
937 ideal h2=id_complement(h), aset=idInit(1,1);
938 int i,deg1=pTotaldegree(S);
939 int tdeg=deg1+ddeg;
940 if(tdeg!=0)
941 {
942 std::vector<int> v,bv=support1(S),in;
943 std::vector<std::vector<int> > hvs=supports(h);
944 ideal ia=id_MaxIdeal(tdeg, currRing);
945 for(i=0;i<IDELEMS(ia);i++)
946 {
947 v=support1(ia->m[i]);
948 in=vecIntersection(v,bv);
949 if(vInvsl(v,hvs)&&in.size()==0)
950 {
951 idInsertPoly(aset, ia->m[i]);
952 }
953 }
954 idSkipZeroes(aset);
955 }
956 else idInsertPoly(aset,e);
957 return(aset);
958}
959
960//returns true if support(p) union support(a) minus support(b) is face,
961//otherwise returns false
962//(the vector version of mabcondition)
963static bool mabconditionv(std::vector<std::vector<int> > hvs,std::vector<int> pv,std::vector<int> av,std::vector<int> bv)
964{
965 std::vector<int> uv=vecUnion(pv,av);
966 uv=vecMinus(uv,bv);
967 if(vInvsl(uv,hvs))
968 {
969 return(true);
970 }
971 return(false);
972}
973
974// returns the set of nonfaces p where mabconditionv(h, p, a, b) is true
975static std::vector<std::vector<int> > Mabv(ideal h,poly a,poly b)
976{
977 std::vector<int> av=support1(a), bv=support1(b), pv, vec;
978 ideal h2=id_complement(h);
979 std::vector<std::vector<int> > hvs=supports(h), h2v=supports(h2), vecs;
980 for(unsigned i=0;i<h2v.size();i++)
981 {
982 pv=h2v[i];
983 if(mabconditionv(hvs,pv,av,bv))
984 {
985 vecs.push_back(pv);
986 }
987 }
988 return vecs;
989}
990
991/***************************************************************************/
992//For solving the equations which has form of x_i-x_j.(equations got from T_1)
993/***************************************************************************/
994
995//subroutine for soleli1
996static std::vector<int> eli1(std::vector<int> eq1,std::vector<int> eq2)
997{
998 int i,j;
999 std::vector<int> eq;
1000 if(eq1[0]==eq2[0])
1001 {
1002 i=eq1[1];j=eq2[1];
1003 eq.push_back(i);
1004 eq.push_back(j);
1005 }
1006 else
1007 {
1008 eq=eq2;
1009 }
1010 return(eq);
1011}
1012
1013/*
1014//get triangular form(eqs.size()>0)
1015static std::vector<std::vector<int> > soleli1( std::vector<std::vector<int> > eqs)
1016{
1017 int i,j;
1018 std::vector<int> yaya;
1019 std::vector<std::vector<int> > pre=eqs, ppre, re;
1020 if(eqs.size()>0)
1021 {
1022 re.push_back(eqs[0]);
1023 pre.erase(pre.begin());
1024 }
1025 for(i=0;i<re.size(),pre.size()>0;i++)
1026 {
1027 yaya=eli1(re[i],pre[0]);
1028 re.push_back(yaya);
1029 for(j=1;j<pre.size();j++)
1030 {
1031 ppre.push_back(eli1(re[i],pre[j]));
1032 }
1033 pre=ppre;
1034 ppre.resize(0);
1035 }
1036 return re;
1037}*/
1038//make sure the first element is smaller that the second one
1039static std::vector<int> keeporder( std::vector<int> vec)
1040{
1041 std::vector<int> yaya;
1042 int n;
1043 if(vec[0]>vec[1])
1044 {
1045 n=vec[0];
1046 vec[0]=vec[1];
1047 vec[1]=n;
1048 }
1049 return vec;
1050}
1051
1052static std::vector<std::vector<int> > soleli1( std::vector<std::vector<int> > eqs)
1053{
1054 int i;
1055 std::vector<int> yaya;
1056 std::vector<std::vector<int> > pre=eqs, ppre, re;
1057 if(eqs.size()>0)
1058 {
1059 re.push_back(eqs[0]);
1060 pre.erase(pre.begin());
1061 }
1062 while(pre.size()>0)
1063 {
1064 yaya=keeporder(eli1(re[0],pre[0]));
1065 for(i=1;i<re.size();i++)
1066 {
1067 if(!vInvsl(yaya, re))
1068 {
1069 yaya=eli1(re[i],yaya);
1070 yaya=keeporder(yaya);
1071 }
1072 }
1073 if(!vInvsl(yaya, re))
1074 {
1075 re.push_back(yaya);
1076 }
1077 pre.erase(pre.begin());
1078 }
1079 return re;
1080}
1081
1082// input is a set of equations who is of triangular form(every equations has a form of x_i-x_j)
1083// n is the number of variables
1084//get the free variables and the dimension
1085static std::vector<int> freevars(int n, std::vector<int> bset, std::vector<std::vector<int> > gset)
1086{
1087 int ql=gset.size(), bl=bset.size(), i;
1088 std::vector<int> mvar, fvar;
1089 for(i=0;i<bl;i++)
1090 {
1091 mvar.push_back(bset[i]);
1092 }
1093 for(i=0;i<ql;i++)
1094 {
1095 mvar.push_back(gset[i][0]);
1096 }
1097 for(i=1;i<=n;i++)
1098 {
1099 if(!IsinL(i,mvar))
1100 {
1101 fvar.push_back(i);
1102 }
1103 }
1104 return fvar;
1105}
1106
1107//return the set of free variables except the vnum one
1108static std::vector<int> fvarsvalue(int vnum, std::vector<int> fvars)
1109{
1110 int i;
1111 std::vector<int> fset=fvars;
1112 for(i=0;i<fset.size();i++)
1113 {
1114 if(fset[i]==vnum)
1115 {
1116 fset.erase(fset.begin()+i);
1117 break;
1118 }
1119 }
1120 return fset;
1121}
1122
1123//returns the simplified bset and gset
1124//enlarge bset, simplify gset
1125static std::vector<std::vector<int> > vAbsorb( std::vector<int> bset,std::vector<std::vector<int> > gset)
1126{
1127 std::vector<int> badset=bset;
1128 int i,j,m, bl=bset.size(), gl=gset.size();
1129 for(i=0;i<bl;i++)
1130 {
1131 m=badset[i];
1132 for(j=0;j<gl;j++)
1133 {
1134 if(gset[j][0]==m && !IsinL(gset[j][1],badset))
1135 {
1136 badset.push_back(gset[j][1]);
1137 gset.erase(gset.begin()+j);
1138 j--;
1139 gl--;
1140 bl++;
1141 }
1142 else if(!IsinL(gset[j][0],badset) && gset[j][1]==m)
1143 {
1144 badset.push_back(gset[j][0]);
1145 gset.erase(gset.begin()+j);
1146 j--;
1147 gl--;
1148 bl++;
1149 }
1150 else if(IsinL(gset[j][0],badset) && IsinL(gset[j][1],badset))
1151 {
1152 gset.erase(gset.begin()+j);
1153 j--;
1154 gl--;
1155 }
1156 else
1157 {
1158 ;
1159 }
1160 }
1161 }
1162 if(badset.size()==0) badset.push_back(0);
1163 gset.push_back(badset);
1164 return gset;
1165}
1166
1167//the labels of new variables are started with 1
1168//returns a vector of solution space according to index
1169static std::vector<int> vecbase1(int num, std::vector<int> oset)
1170{
1171 int i;
1172 std::vector<int> base;
1173 for(i=0;i<num;i++)
1174 {
1175 if(IsinL(i+1,oset))
1176 base.push_back(1);
1177 else
1178 base.push_back(0);
1179 }
1180 return base;
1181}
1182
1183//returns a vector which has length of n,
1184//and all the entries are 0.
1185static std::vector<int> make0(int n)
1186{
1187 int i;
1188 std::vector<int> vec;
1189 for(i=0;i<n;i++)
1190 {
1191 vec.push_back(0);
1192 }
1193 return vec;
1194}
1195
1196//returns a vector which has length of n,
1197//and all the entries are 1.
1198static std::vector<int> make1(int n)
1199{
1200 int i;
1201 std::vector<int> vec;
1202 for(i=0;i<n;i++)
1203 {
1204 vec.push_back(1);
1205 }
1206 return vec;
1207}
1208
1209//input gset must be the triangular form after zero absorbing according to the badset,
1210//bset must be the zero set after absorbing.
1211static std::vector<int> ofindbases1(int num, int vnum, std::vector<int> bset,std::vector<std::vector<int> > gset)
1212{
1213 std::vector<std::vector<int> > goodset;
1214 std::vector<int> fvars=freevars(num, bset, gset), oset, base;
1215 std::vector<int> zset=fvarsvalue(vnum, fvars);
1216 zset=vecUnion(zset,bset);
1217 oset.push_back(vnum);
1218 goodset=vAbsorb(oset, gset);
1219 oset=goodset[goodset.size()-1];
1220 goodset.erase(goodset.end());
1221 base= vecbase1(num, oset);
1222 return base;
1223}
1224
1225//input gset must be the triangular form after zero absorbing according to the badset
1226//bset must be the zero set after absorbing
1227static std::vector<std::vector<int> > ofindbases(int num, std::vector<int> bset,std::vector<std::vector<int> > gset)
1228{
1229 int i,m;
1230 std::vector<std::vector<int> > bases;
1231 std::vector<int> fvars=freevars(num, bset, gset), base1;
1232 if (fvars.size()==0)
1233 {
1234 base1=make0(num);
1235 bases.push_back(base1);
1236 }
1237 else
1238 {
1239 for(i=0;i<fvars.size();i++)
1240 {
1241 m=fvars[i];
1242 base1=ofindbases1(num, m, bset, gset);
1243 bases.push_back(base1);
1244 }
1245 }
1246 //PrintS("They are the bases for the solution space:\n");
1247 //listsprint(bases);
1248 return bases;
1249}
1250
1251//gset is a set of equations which have forms of x_i-x_j
1252//num is the number of varialbes also the length of the set which we need to consider
1253//output is trigular form of gset and badset where x_i=0
1254static std::vector<std::vector<int> > eli2(int num, std::vector<int> bset,std::vector<std::vector<int> > gset)
1255{
1256 std::vector<int> badset;
1257 std::vector<std::vector<int> > goodset, solve;
1258//PrintS("This is the input bset\n");listprint(bset);
1259//PrintS("This is the input gset\n");listsprint(gset);
1260 if(gset.size()!=0)//gset is not empty
1261 {
1262 //find all the variables which are zeroes
1263
1264 if(bset.size()!=0)//bset is not empty
1265 {
1266 goodset=vAbsorb(bset, gset);//e.g. x_1=0, put x_i into the badset if x_i-x_1=0 or x_1-x_i=0
1267 int m=goodset.size();
1268 badset=goodset[m-1];
1269 goodset.erase(goodset.end());
1270 }
1271 else //bset is empty
1272 {
1273 goodset=gset;//badset is empty
1274 }//goodset is already the set which doesn't contain zero variables
1275//PrintS("This is the badset after absorb \n");listprint(badset);
1276//PrintS("This is the goodset after absorb \n");listsprint(goodset);
1277 goodset=soleli1(goodset);//get the triangular form of goodset
1278//PrintS("This is the goodset after triangulization \n");listsprint(goodset);
1279 solve=ofindbases(num,badset,goodset);
1280 }
1281 else
1282 {
1283 solve=ofindbases(num,bset,gset);
1284 }
1285//PrintS("This is the solution\n");listsprint(solve);
1286 return solve;
1287}
1288
1289/********************************************************************/
1290/************************links***********************************/
1291
1292//returns the links of face a in simplicial complex X
1293static std::vector<std::vector<int> > links(poly a, ideal h)
1294{
1295 int i;
1296 std::vector<std::vector<int> > lk,X=supports(h);
1297 std::vector<int> U,In,av=support1(a);
1298 for(i=0;i<X.size();i++)
1299 {
1300 U=vecUnion(av,X[i]);
1301 In=vecIntersection(av,X[i]);
1302 if( In.size()==0 && vInvsl(U,X))
1303 {
1304 //PrintS("The union of them is FACE and intersection is EMPTY!\n");
1305 lk.push_back(X[i]);
1306 }
1307 else
1308 {
1309 ;
1310 }
1311 }
1312 return lk;
1313}
1314
1315static int redefinedeg(poly p, int num)
1316{
1317 int deg=0, deg0;
1318 for(int i=1;i<=currRing->N;i++)
1319 {
1320 deg0=pGetExp(p, i);
1321 if(i>num)
1322 {
1323 deg= deg+2*deg0;
1324 }
1325 else
1326 {
1327 deg=deg+deg0;
1328 }
1329 }
1330 //Print("the new degree is: %d\n", deg);
1331 return (deg);
1332}
1333
1334// the degree of variables should be same
1335static ideal p_a(ideal h)
1336{
1337 poly p;
1338 int i,j,deg=0,deg0;
1339 ideal aset=idCopy(h),ia,h1=idsrRing(h);
1340//PrintS("idsrRing is:\n");id_print(h1);
1341 std::vector<int> as;
1342 std::vector<std::vector<int> > hvs=supports(h);
1343 for(i=0;i<IDELEMS(h1);i++)
1344 {
1345 deg0=pTotaldegree(h1->m[i]);
1346 if(deg < deg0)
1347 deg=deg0;
1348 }
1349 for(i=2;i<=deg;i++)
1350 {
1351 ia=id_MaxIdeal(i, currRing);
1352 for(j=0;j<IDELEMS(ia);j++)
1353 {
1354 p=pCopy(ia->m[j]);
1355 if(!IsInX(p,h))
1356 {
1357 as=support1(p);
1358 if(vInvsl(as,hvs))
1359 {
1360 idInsertPoly(aset, p);
1361 }
1362 }
1363 }
1364 }
1365 idSkipZeroes(aset);
1366 return(aset);
1367}
1368
1369/*only for the exampels whose variables has degree more than 1*/
1370/*ideal p_a(ideal h)
1371{
1372 poly e=pOne(), p;
1373 int i,j,deg=0,deg0, ord=4;
1374 ideal aset=idCopy(h),ia,h1=idsrRing(h);
1375//PrintS("idsrRing is:\n");id_print(h1);
1376 std::vector<int> as;
1377 std::vector<std::vector<int> > hvs=supports(h);
1378 for(i=0;i<IDELEMS(h1);i++)
1379 {
1380 deg0=redefinedeg(h1->m[i],ord);
1381 if(deg < deg0)
1382 deg=deg0;
1383 }
1384 for(i=2;i<=deg;i++)
1385 {
1386 ia=id_MaxIdeal(i, currRing);
1387 for(j=0;j<IDELEMS(ia);j++)
1388 {
1389 p=pCopy(ia->m[j]);
1390 if(!IsInX(p,h))
1391 {
1392 as=support1(p);
1393 if(vInvsl(as,hvs))
1394 {
1395 idInsertPoly(aset, p);
1396 }
1397 }
1398 }
1399 }
1400 idSkipZeroes(aset);
1401 return(aset);
1402}*/
1403
1404static std::vector<int> vertset(std::vector<std::vector<int> > vecs)
1405{
1406 int i,j;
1407 std::vector<int> vert;
1408 std::vector<std::vector<int> > vvs;
1409 for(i=1;i<=currRing->N;i++)
1410 {
1411 for(j=0;j<vecs.size();j++)
1412 {
1413 if(IsinL(i, vecs[j]))
1414 {
1415 if(!IsinL(i , vert))
1416 {
1417 vert.push_back(i);
1418 }
1419 break;
1420 }
1421 }
1422 }
1423 return (vert);
1424}
1425
1426//smarter way
1427static ideal p_b(ideal h, poly a)
1428{
1429 std::vector<std::vector<int> > pbv,lk=links(a,h), res;
1430 std::vector<int> vert=vertset(lk), bv;
1431 res=b_subsets(vert);
1432 int i, adg=pTotaldegree(a);
1433 poly e=pOne();
1434 ideal idd=idInit(1,1);
1435 for(i=0;i<res.size();i++)
1436 {
1437 if(res[i].size()==adg)
1438 pbv.push_back(res[i]);
1439 }
1440 if(pEqualPolys(a,e))
1441 {
1442 idInsertPoly(idd, e);
1443 idSkipZeroes(idd);
1444 return (idd);
1445 }
1446 idd=idMaken(pbv);
1447 return(idd);
1448}
1449
1450/*//dump way to get pb
1451// the degree of variables should be same
1452static ideal p_b(ideal h, poly a)
1453{
1454 std::vector<std::vector<int> > pbv,lk=links(a,h),res;
1455// PrintS("Its links are :\n");id_print(idMaken(lk));
1456 res=id_subsets(lk);
1457 //PrintS("res is :\n");listsprint(res);
1458 std::vector<int> bv;
1459 ideal bset=findb(h);
1460 int i,j,nu=res.size(),adg=pTotaldegree(a);
1461 poly e=pOne();ideal idd=idInit(1,1);
1462 for(i=0;i<res.size();i++)
1463 {
1464 if(res[i].size()==adg)
1465 pbv.push_back(res[i]);
1466 }
1467 if(pEqualPolys(a,e)){idInsertPoly(idd, e); idSkipZeroes(idd); return (idd);}
1468 for(i=0;i<nu;i++)
1469 {
1470 for(j=i+1;j<nu;j++)
1471 {
1472 if(res[i].size()!=0 && res[j].size()!=0)
1473 {
1474 bv = vecUnion(res[i], res[j]);
1475 if(IsInX(pMaken(bv),bset) && bv.size()==adg && !vInvsl(bv,pbv))
1476 {pbv.push_back(bv);}
1477 }
1478 }
1479 }
1480 idd=idMaken(pbv);
1481 //id_print(idd);
1482 return(idd);
1483}*/
1484
1485// also only for the examples whose variables have degree more than 1(ndegreeb and p_b)
1486/*int ndegreeb(std::vector<int> vec, int num)
1487{
1488 int deg, deg0=0;
1489 for(int i=0;i<vec.size();i++)
1490 {
1491 if(vec[i]>num)
1492 {
1493 deg0++;
1494 }
1495 }
1496 deg=vec.size()+deg0;
1497 return(deg);
1498}
1499
1500static ideal p_b(ideal h, poly a)
1501{
1502 std::vector<std::vector<int> > pbv,lk=links(a,h),res;
1503// PrintS("Its links are :\n");id_print(idMaken(lk));
1504 res=id_subsets(lk);
1505 //PrintS("res is :\n");listsprint(res);
1506 std::vector<int> bv;
1507 ideal bset=findb(h);
1508 int i,j,nu=res.size(),ord=4,adg=redefinedeg(a, ord);
1509 poly e=pOne();ideal idd=idInit(1,1);
1510 for(i=0;i<res.size();i++)
1511 {
1512 if(ndegreeb(res[i],ord)==adg)
1513 pbv.push_back(res[i]);
1514 }
1515 if(pEqualPolys(a,e)){idInsertPoly(idd, e); idSkipZeroes(idd); return (idd);}
1516 for(i=0;i<nu;i++)
1517 {
1518 for(j=i+1;j<nu;j++)
1519 {
1520 if(res[i].size()!=0 && res[j].size()!=0)
1521 {
1522 bv = vecUnion(res[i], res[j]);
1523 //PrintS("bv is :\n");listprint(bv);
1524 //Print("bv's degree is : %d\n", ndegreeb(bv,ord));
1525 if(IsInX(pMaken(bv),bset) && ndegreeb(bv,ord)==adg && !vInvsl(bv,pbv))
1526 {
1527 pbv.push_back(bv);
1528 }
1529 }
1530 }
1531 }
1532 idd=idMaken(pbv);
1533 //id_print(idd);
1534 return(idd);
1535}*/
1536
1537//input is a squarefree monomial p
1538//output is all the squarefree monomials which could divid p(including p itself?)
1539static ideal psubset(poly p)
1540{
1541 int i,j,max=pTotaldegree(p);
1542 ideal h1,mons, id_re=idInit(1,1);
1543 for(i=1;i<max;i++)
1544 {
1545 mons=id_MaxIdeal(i, currRing);
1546 h1=sfreemon(mons,i);
1547 for(j=0;j<IDELEMS(h1);j++)
1548 {
1549 if(p_DivisibleBy(h1->m[j],p,currRing))
1550 idInsertPoly(id_re, h1->m[j]);
1551 }
1552 }
1553 idSkipZeroes(id_re);
1554 //PrintS("This is the facset\n");
1555 //id_print(id_re);
1556 return id_re;
1557}
1558
1559//inserts a new vector which has two elements a and b into vector gset (which is a vector of vectors)
1560//(especially for gradedpiece1 and gradedpiece1n)
1561static std::vector<std::vector<int> > listsinsertlist(std::vector<std::vector<int> > gset, int a, int b)
1562{
1563 std::vector<int> eq;
1564 eq.push_back(a);
1565 eq.push_back(b);
1566 gset.push_back(eq);
1567 return gset;
1568}
1569
1570static std::vector<int> makeequation(int i,int j, int t)
1571{
1572 std::vector<int> equation;
1573 equation.push_back(i);
1574 equation.push_back(j);
1575 equation.push_back(t);
1576 //listprint(equation);
1577 return equation;
1578}
1579
1580/****************************************************************/
1581//only for solving the equations obtained from T^2
1582//input should be a vector which has only 3 entries
1583static poly pMake3(std::vector<int> vbase)
1584{
1585 int co=1;
1586 poly p,q=0;
1587 for(int i=0;i<3;i++)
1588 {
1589 if(vbase[i]!=0)
1590 {
1591 if(i==1) co=-1;
1592 p = pOne();pSetExp(p, vbase[i], 1);pSetm(p);pSetCoeff(p, nInit(co));
1593 }
1594 else p=0;
1595 q = pAdd(q, p);
1596 co=1;
1597 }
1598 return q;
1599}
1600
1601static ideal idMake3(std::vector<std::vector<int> > vecs)
1602{
1603 ideal id_re=idInit(1,1);
1604 poly p;
1605 int i,lv=vecs.size();
1606 for(i=0;i<lv;i++)
1607 {
1608 p=pMake3(vecs[i]);
1609 idInsertPoly(id_re, p);
1610 }
1611 idSkipZeroes(id_re);
1612 return id_re;
1613}
1614
1615/****************************************************************/
1616
1617//change the current ring to a new ring which is in num new variables
1618static void equmab(int num)
1619{
1620 int i;
1621 //Print("There are %d new variables for equations solving.\n",num);
1622 ring r=currRing;
1623 char** tt;
1624 coeffs cf=nCopyCoeff(r->cf);
1625 tt=(char**)omAlloc(num*sizeof(char *));
1626 for(i=0; i <num; i++)
1627 {
1628 tt[i] = (char*)omalloc(10); //if required enlarge it later
1629 sprintf (tt[i], "t(%d)", i+1);
1630 tt[i]=omStrDup(tt[i]);
1631 }
1632 ring R=rDefault(cf,num,tt,ringorder_lp);
1634 IDRING(h)=rCopy(R);
1635 rSetHdl(h);
1636}
1637
1638//returns the trivial case of T^1
1639//b must only contain one variable
1640static std::vector<int> subspace1(std::vector<std::vector<int> > mv, std::vector<int> bv)
1641{
1642 int i, num=mv.size();
1643 std::vector<int> base;
1644 for(i=0;i<num;i++)
1645 {
1646 if(IsinL(bv[0],mv[i]))
1647 base.push_back(1);
1648 else
1649 base.push_back(0);
1650 }
1651 return base;
1652}
1653
1654/***************************only for T^2*************************************/
1655//vbase only has two elements which records the position of the monomials in mv
1656
1657static std::vector<poly> pMakei(std::vector<std::vector<int> > mv,std::vector<int> vbase)
1658{
1659 poly p;
1660 std::vector<poly> h1;
1661 int n=vbase.size();
1662 for(int i=0;i<n;i++)
1663 {
1664 p=pMaken(mv[vbase[i]]);
1665 h1.push_back(p);
1666 }
1667 return h1;
1668}
1669
1670// returns a ideal according to a set of supports
1671static std::vector<std::vector<poly> > idMakei(std::vector<std::vector<int> > mv,std::vector<std::vector<int> > vecs)
1672{
1673 int i,lv=vecs.size();
1674 std::vector<std::vector<poly> > re;
1675 std::vector<poly> h;
1676 for(i=0;i<lv;i++)
1677 {
1678 h=pMakei(mv,vecs[i]);
1679 re.push_back(h);
1680 }
1681 //PrintS("This is the metrix M:\n");
1682 //listsprint(vecs);
1683 //PrintS("the ideal according to metrix M is:\n");
1684 return re;
1685}
1686
1687/****************************************************************/
1688
1689//return the graded pieces of cohomology T^1 according to a,b
1690//original method (only for debugging)
1691static void gradedpiece1(ideal h,poly a,poly b)
1692{
1693 int i,j,m;
1694 ideal sub=psubset(b);
1695 std::vector<int> av=support1(a), bv=support1(b), bad, vv;
1696 std::vector<std::vector<int> > hvs=supports(h), sbv=supports(sub), mv=Mabv(h,a,b),good;
1697 m=mv.size();
1698 ring r=currRing;
1699 if( m > 0 )
1700 {
1701 for(i=0;i<m;i++)
1702 {
1703 if(!vsubset(bv,mv[i]))
1704 {
1705 bad.push_back(i+1);
1706 }
1707 }
1708 for(i=0;i<m;i++)
1709 {
1710 for(j=i+1;j<m;j++)
1711 {
1712 vv=vecUnion(mv[i],mv[j]);
1713 if(mabconditionv(hvs,vv,av,bv))
1714 {
1715 good=listsinsertlist(good,i+1,j+1);
1716 }
1717 else
1718 {
1719 //PrintS("They are not in Mabt!\n");
1720 ;
1721 }
1722 }
1723 }
1724 std::vector<std::vector<int> > solve=eli2(m,bad,good);
1725 if(bv.size()!=1)
1726 {
1727 //PrintS("This is the solution of coefficients:\n");
1729 }
1730 else
1731 {
1732 std::vector<int> su=subspace1(mv,bv);
1733 //PrintS("This is the solution of subspace:\n");
1734 //listprint(su);
1735 std::vector<std::vector<int> > suu;
1736 suu.push_back(su);
1737 equmab(solve[0].size());
1738 std::vector<std::vector<int> > solves=vecqring(solve,suu);
1739 //PrintS("This is the solution of coefficients:\n");
1740 listsprint(solves);
1741 rChangeCurrRing(r);
1742 }
1743 }
1744 else
1745 {
1746 PrintS("No element considered!\n");
1747 }
1748}
1749
1750//Returns true if b can divide p*q
1751static bool condition1for2(std::vector<int > pv,std::vector<int > qv,std::vector<int > bv)
1752{
1753 std::vector<int > vec=vecUnion(pv,qv);
1754 if(vsubset(bv,vec))
1755 {
1756 //PrintS("condition1for2 yes\n");
1757 return true;
1758 }
1759 //PrintS("condition1for2 no\n");
1760 return false;
1761}
1762
1763//Returns true if support(p) union support(q) union support(s) union support(a) minus support(b) is face
1764static bool condition2for2(std::vector<std::vector<int> > hvs, std::vector<int> pv, std::vector<int> qv, std::vector<int> sv, std::vector<int> av, std::vector<int> bv)
1765{
1766 std::vector<int> vec=vecUnion(pv,qv);
1767 vec=vecUnion(vec,sv);
1768 if(mabconditionv(hvs,vec,av,bv))
1769 {
1770 //PrintS("condition2for2 yes\n");
1771 return (true);
1772 }
1773 //PrintS("condition2for2 no\n");
1774 return (false);
1775}
1776
1777static bool condition3for2(std::vector<std::vector<int> > hvs, std::vector<int> pv, std::vector<int> qv, std::vector<int> av, std::vector<int> bv)
1778{
1779 std::vector<int> v1,v2,v3;
1780 v1=vecIntersection(pv,qv);//intersection
1781 v2=vecUnion(pv,qv);
1782 v2=vecUnion(v2,av);
1783 v2=vecMinus(v2,bv);
1784 v3=vecUnion(v1,v2);
1785 if(vInvsl(v3,hvs))
1786 {
1787 //PrintS("condition3for2 yes\n");
1788 return(true);
1789 }
1790 //PrintS("condition3for2 no\n");
1791 return(false);
1792}
1793
1794/****************solve the equations got from T^2*********************/
1795
1796static ideal getpresolve(ideal h)
1797{
1798 //ring r=currRing;
1799 //assume (LIB "presolve.lib");
1800 sleftv a;a.Init();
1801 a.rtyp=IDEAL_CMD;a.data=(void*)h;
1802 idhdl solve=ggetid("elimlinearpart");
1803 if(solve==NULL)
1804 {
1805 WerrorS("presolve.lib are not loaded!");
1806 return NULL;
1807 }
1809 //PrintS("no errors here\n");
1810 if(sl)
1811 {
1812 WerrorS("error in solve!");
1813 }
1815 ideal re=(ideal)L->m[4].CopyD();
1816 //iiRETURNEXPR.CleanUp();
1818 //PrintS("no errors here\n");
1819 //idSkipZeroes(re);
1820 //id_print(re);
1821 return re;
1822}
1823
1824static std::vector<int> numfree(ideal h)
1825{
1826 int i,j;
1827 std::vector<int> fvar;
1828 for(j=1;j<=currRing->N;j++)
1829 {
1830 for(i=0;i<IDELEMS(h);i++)
1831 {
1832 if(vInp(j,h->m[i]))
1833 {
1834 fvar.push_back(j);
1835 break;
1836 }
1837 }
1838 }
1839 //Print("There are %d free variables in total\n",num);
1840 return fvar;
1841}
1842
1843static std::vector<std::vector<int> > canonicalbase(int n)
1844{
1845 std::vector<std::vector<int> > vecs;
1846 std::vector<int> vec;
1847 int i,j;
1848 for(i=0;i<n;i++)
1849 {
1850 for(j=0;j<n;j++)
1851 {
1852 if(i==j)
1853 vec.push_back(1);
1854 else
1855 vec.push_back(0);
1856 }
1857 vecs.push_back(vec);
1858 vec.clear();
1859 }
1860 return vecs;
1861}
1862
1863static std::vector<std::vector<int> > getvector(ideal h,int n)
1864{
1865 std::vector<int> vec;
1866 std::vector<std::vector<int> > vecs;
1867 ideal h2=idCopy(h);
1868 if(!idIs0(h))
1869 {
1870 ideal h1=getpresolve(h2);
1871 poly q,e=pOne();
1872 int lg=IDELEMS(h1),n,i,j,t;
1873 std::vector<int> fvar=numfree(h1);
1874 n=fvar.size();
1875 if(n==0)
1876 {
1877 vec=make0(IDELEMS(h1));vecs.push_back(vec);//listsprint(vecs);
1878 }
1879 else
1880 {
1881 for(t=0;t<n;t++)
1882 {
1883 vec.clear();
1884 for(i=0;i<lg;i++)
1885 {
1886 q=pCopy(h1->m[i]);
1887 //pWrite(q);
1888 if(q==0)
1889 {
1890 vec.push_back(0);
1891 }
1892 else
1893 {
1894 q=p_Subst(q, fvar[t], e,currRing);
1895 //Print("the %dth variable was substituted by 1:\n",fvar[t]);
1896 //pWrite(q);
1897 for(j=0;j<n;j++)
1898 {
1899 //Print("the %dth variable was substituted by 0:\n",fvar[j]);
1900 q=p_Subst(q, fvar[j],0,currRing);
1901 //pWrite(q);
1902 }
1903 if(q==0)
1904 {
1905 vec.push_back(0);
1906 }
1907 else
1908 {
1909 vec.push_back(n_Int(pGetCoeff(q),currRing->cf));
1910 }
1911 }
1912 }
1913 //listprint(vec);
1914 vecs.push_back(vec);
1915 }
1916 }
1917 }
1918 else
1919 {vecs=canonicalbase(n);}
1920 //listsprint(vecs);
1921 return vecs;
1922}
1923
1924/**************************************************************************/
1925
1926//subspace of T2(find all the possible values of alpha)
1927static std::vector<int> findalpha(std::vector<std::vector<int> > mv, std::vector<int> bv)
1928{
1929 std::vector<int> alset;
1930 for(unsigned i=0;i<mv.size();i++)
1931 {
1932 if(vsubset(bv,mv[i]))
1933 {
1934 alset.push_back(i);
1935 }
1936 }
1937 //Print("This is the alpha set, and the subspace is dim-%ld\n",alset.size());
1938 //listprint(alset);
1939 return alset;
1940}
1941
1942static std::vector<int> subspacet1(int num, std::vector<std::vector<int> > ntvs)
1943{
1944 int i, j, t, n=ntvs.size();
1945 std::vector<int> subase;
1946 for(t=0;t<n;t++)
1947 {
1948 i=ntvs[t][0];
1949 j=ntvs[t][1];
1950 if(i==(num))
1951 {
1952 subase.push_back(1);
1953 }
1954 else if(j==num)
1955 {
1956 subase.push_back(-1);
1957 }
1958 else
1959 {
1960 subase.push_back(0);
1961 }
1962 }
1963 //Print("This is the basis w.r.t. %dth polynomial in alpha set\n",num);
1964 //listprint(subase);
1965 return subase;
1966}
1967
1968//subspace for T^2(mab method)
1969static std::vector<std::vector<int> > subspacet(std::vector<std::vector<int> > mv, std::vector<int> bv,std::vector<std::vector<int> > ntvs)
1970{
1971 std::vector<int> alset=findalpha(mv,bv), subase;
1972 std::vector<std::vector<int> > subases;
1973 for(unsigned i=0;i<alset.size();i++)
1974 {
1975 subase=subspacet1(alset[i],ntvs);
1976 subases.push_back(subase);
1977 }
1978 //PrintS("These are the bases for the subspace:\n");
1979 //listsprint(subases);
1980 return subases;
1981}
1982
1983static std::vector<std::vector<int> > mabtv(std::vector<std::vector<int> > hvs, std::vector<std::vector<int> > Mv, std::vector<int> av, std::vector<int> bv)
1984{
1985 std::vector<int> v1,var;
1986 std::vector<std::vector<int> > vars;
1987 for(unsigned i=0;i<Mv.size();i++)
1988 {
1989 for(unsigned j=i+1;j<Mv.size();j++)
1990 {
1991 var.clear();
1992 v1=vecUnion(Mv[i],Mv[j]);
1993 if(mabconditionv(hvs, v1, av, bv))
1994 {
1995 var.push_back(i);
1996 var.push_back(j);
1997 vars.push_back(var);
1998 }
1999 }
2000 }
2001 return vars;
2002}
2003
2004//fix the problem of the number of the new variables
2005//original method for T^2(only for debugging)
2006static void gradedpiece2(ideal h,poly a,poly b)
2007{
2008 int t0,t1,t2,i,j,t,m;
2009 ideal sub=psubset(b);
2010 ring r=rCopy(currRing);
2011 std::vector<std::vector<int> > hvs=supports(h), mv=Mabv(h,a,b), mts, vecs,vars;
2012 std::vector<int> av=support1(a), bv=support1(b), vec,var;
2013 mts=mabtv(hvs,mv,av,bv);
2014 PrintS("The homomorphism should map onto:\n");
2015 lpsprint(idMakei(mv,mts));
2016 m=mv.size();
2017 if(m > 0)
2018 {
2019 vars=mabtv(hvs,mv,av,bv);
2020 int vn=vars.size();
2021 for(t0=0;t0<vars.size();t0++)
2022 {
2023 i=vars[t0][0];
2024 j=vars[t0][1];
2025 if(!condition1for2(mv[i],mv[j],bv))//condition 1
2026 {
2027 //PrintS("And they satisfy the condition 1.\n");
2028 vec=makeequation(t0+1,0,0);
2029 //PrintS("So the equation:\n");
2030 //pWrite(p);
2031 //PrintS("holds.\n");
2032 vecs.push_back(vec);
2033 vec.clear();
2034 }
2035 if(condition3for2(hvs,mv[i],mv[j],av,bv))//condition 3
2036 {
2037 //PrintS("And they satisfy the condition 3.\n");
2038 vec=makeequation(t0+1,0,0);
2039 //PrintS("So the equation: \n");
2040 //pWrite(p);
2041 //PrintS("holds.\n");
2042 vecs.push_back(vec);
2043 vec.clear();
2044 }
2045 for(t1=t0+1;t1<vars.size();t1++)
2046 {
2047 for(t2=t1+1;t2<vars.size();t2++)
2048 {
2049 if(vars[t0][0]==vars[t1][0]&&vars[t1][1]==vars[t2][1]&&vars[t0][1]==vars[t2][0])
2050 {
2051 i=vars[t0][0];
2052 j=vars[t0][1];
2053 t=vars[t1][1];
2054 if(condition2for2(hvs,mv[i],mv[j],mv[t],av,bv))//condition 2
2055 {
2056 vec=makeequation(t0+1,t1+1,t2+1);
2057 vecs.push_back(vec);
2058 vec.clear();
2059 }
2060 }
2061 }
2062 }
2063 }
2064 //PrintS("this is EQUATIONS:\n");
2065 //listsprint(vecs);
2066 equmab(vn);
2067 ideal id_re=idMake3(vecs);
2068 //id_print(id_re);
2069 std::vector<std::vector<int> > re=getvector(id_re,vn);
2070 PrintS("this is the solution for ideal :\n");
2071 listsprint(re);
2072 rChangeCurrRing(r);
2073 std::vector<std::vector<int> > sub=subspacet(mv, bv,vars);
2074 PrintS("this is the solution for subspace:\n");
2075 listsprint(sub);
2076 equmab(vn);
2077 std::vector<std::vector<int> > solve=vecqring(re, sub);
2078 PrintS("This is the solution of coefficients:\n");
2080 rChangeCurrRing(r);
2081 }
2082 else
2083 {
2084 PrintS("No element considered!");
2085 }
2086}
2087
2088/**********************************************************************/
2089//For the method of N_{a-b}
2090
2091//returns true if pv(support of monomial) satisfies pv union av minus bv is in hvs
2092static bool nabconditionv(std::vector<std::vector<int> > hvs, std::vector<int> pv, std::vector<int> av, std::vector<int> bv)
2093{
2094 std::vector<int> vec1=vecIntersection(pv,bv), vec2=vecUnion(pv,bv);
2095 int s1=vec1.size();
2096 if(!vInvsl(vec2,hvs) && s1==0 && vsubset(av,pv))
2097 {
2098 //PrintS("nab condition satisfied\n");
2099 return(true);
2100 }
2101 //PrintS("nab condition not satisfied\n");
2102 return(false);
2103}
2104
2105//returns N_{a-b}
2106static std::vector<std::vector<int> > Nabv(std::vector<std::vector<int> > hvs, std::vector<int> av, std::vector<int> bv)
2107{
2108 std::vector<std::vector<int> > vecs;
2109 int num=hvs.size();
2110 for(int i=0;i<num;i++)
2111 {
2112 if(nabconditionv(hvs,hvs[i],av,bv))
2113 {
2114 //PrintS("satisfy:\n");
2115 vecs.push_back(hvs[i]);
2116 }
2117 }
2118 return vecs;
2119}
2120
2121//returns true if pv union qv union av minus bv is in hvs
2122//hvs is simplicial complex
2123static bool nabtconditionv(std::vector<std::vector<int> > hvs, std::vector<int> pv, std::vector<int> qv)
2124{
2125 std::vector<int> v1;
2126 v1=vecUnion(pv,qv);
2127 if(vInvsl(v1,hvs))
2128 {
2129 return (true);
2130 }
2131 return (false);
2132}
2133
2134//returns N_{a-b}^(2)
2135static std::vector<std::vector<int> > nabtv(std::vector<std::vector<int> > hvs, std::vector<std::vector<int> > Nv, std::vector<int> av, std::vector<int> bv)
2136{
2137 std::vector<int> v1,var;
2138 std::vector<std::vector<int> > vars;
2139 for(unsigned i=0;i<Nv.size();i++)
2140 {
2141 for(unsigned j=i+1;j<Nv.size();j++)
2142 {
2143 var.clear();
2144 if(nabtconditionv(hvs, Nv[i], Nv[j]))
2145 {
2146 var.push_back(i);
2147 var.push_back(j);
2148 vars.push_back(var);
2149 }
2150 }
2151 }
2152 return vars;
2153}
2154
2155//p must be the monomial which is a face
2156// ideal sub=psubset(b); bvs=supports(sub);
2157static bool tNab(std::vector<std::vector<int> > hvs, std::vector<int> pv, std::vector<std::vector<int> > bvs)
2158{
2159 std::vector<int> sv;
2160 if(bvs.size()<=1) return false;
2161 for(unsigned i=0;i<bvs.size();i++)
2162 {
2163 sv=vecUnion(pv,bvs[i]);
2164 if(!vInvsl(sv,hvs))
2165 {
2166 return true;
2167 }
2168 }
2169 return false;
2170}
2171
2172static std::vector<int> tnab(std::vector<std::vector<int> > hvs,std::vector<std::vector<int> > nvs,std::vector<std::vector<int> > bvs)
2173{
2174 std::vector<int> pv, vec;
2175 for(unsigned j=0;j<nvs.size();j++)
2176 {
2177 pv=nvs[j];
2178 if(tNab(hvs, pv, bvs))
2179 {
2180 vec.push_back(j);
2181 }
2182 }
2183 return vec;
2184}
2185
2186//the image phi(pv)=pv union av minus bv
2187static std::vector<int> phimage(std::vector<int> pv, std::vector<int> av, std::vector<int> bv)
2188{
2189 std::vector<int> qv=vecUnion(pv,av);
2190 qv=vecMinus(qv,bv);
2191 return qv;
2192}
2193
2194//mvs and nvs are the supports of ideal Mab and Nab
2195//vecs is the solution of nab
2196static std::vector<std::vector<int> > value1(std::vector<std::vector<int> > mvs, std::vector<std::vector<int> > nvs, std::vector<std::vector<int> > vecs,std::vector<int> av, std::vector<int> bv)
2197{
2198 int j;
2199 std::vector<int> pv, base;
2200 std::vector<std::vector<int> > bases;
2201 for(unsigned t=0;t<vecs.size();t++)
2202 {
2203 for(unsigned i=0;i<mvs.size();i++)
2204 {
2205 pv=phimage(mvs[i],av,bv);
2206 for( j=0;j<nvs.size();j++)
2207 {
2208 if(vEvl(pv,nvs[j]))
2209 {
2210 base.push_back(vecs[t][j]);
2211 break;
2212 }
2213 }
2214 if(j==nvs.size())
2215 {
2216 base.push_back(0);
2217 }
2218 }
2219 if(base.size()!=mvs.size())
2220 {
2221 //WerrorS("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1");
2222 WerrorS("Errors in Equations solving (Values Finding)!");
2223 usleep(1000000);
2224 assert(false);
2225
2226 }
2227 bases.push_back(base);
2228 base.clear();
2229 }
2230 return bases;
2231}
2232
2233static intvec *Tmat(std::vector<std::vector<int> > vecs)
2234{
2235 //std::vector<std::vector<int> > solve=gradedpiece1n(h,a,b);
2236 //Print("the size of solve is: %ld\n",solve.size());
2237 //vtm(solve);
2238 intvec *m;
2239 int i,j, a=vecs.size();
2240 if(a==0)
2241 {
2242 m=new intvec(1,1,10);
2243 }
2244 else
2245 {
2246 int b=vecs[0].size();
2247 m=new intvec(a,b,0);
2248 for(i=1;i<=a;i++)
2249 {
2250 for(j=1;j<=b;j++)
2251 {
2252 IMATELEM(*m,i,j)=vecs[i-1][j-1];
2253 }
2254 }
2255 }
2256 return (m);
2257}
2258
2259//returns the set of position number of minimal gens in M
2260static std::vector<int> gensindex(ideal M, ideal ids)
2261{
2262 int i;
2263 std::vector<int> vec,index;
2264 if(!idIs0(M))
2265 {
2266 std::vector<std::vector<int> > vecs=supports(ids);
2267 for(i=0;i<IDELEMS(M);i++)
2268 {
2269 vec=support1(M->m[i]);
2270 if(vInvsl(vec,vecs))
2271 index.push_back(i);
2272 }
2273 }
2274 return (index);
2275}
2276
2277static ideal mingens(ideal h, poly a, poly b)
2278{
2279 int i;
2280 std::vector<std::vector<int> > mv=Mabv(h,a,b);
2281 ideal M=idMaken(mv), hi=idInit(1,1);
2282 std::vector<int> index = gensindex(M, idsrRing(h));
2283 for(i=0;i<index.size();i++)
2284 {
2285 idInsertPoly(hi,M->m[index[i]]);
2286 }
2287 idSkipZeroes(hi);
2288 return (hi);
2289}
2290
2291static std::vector<std::vector<int> > minisolve(std::vector<std::vector<int> > solve, std::vector<int> index)
2292{
2293 int i,j;
2294 std::vector<int> vec,solm;
2295 std::vector<std::vector<int> > solsm;
2296 for(i=0;i<solve.size();i++)
2297 {
2298 vec=solve[i];
2299 for(j=0;j<vec.size();j++)
2300 {
2301 if(IsinL(j,index))
2302 solm.push_back(vec[j]);
2303 }
2304 solsm.push_back(solm);
2305 solm.clear();
2306 }
2307 return (solsm);
2308}
2309
2310//T_1 graded piece(N method)
2311//frame of the most efficient version
2312//regardless of links
2313static intvec * gradedpiece1n(ideal h,poly a,poly b)
2314{
2315 int i,j,co,n;
2316 std::vector<std::vector<int> > hvs=supports(h),mv=Mabv(h,a,b),sbv,nv,good,solve;
2317 std::vector<int> av=support1(a), bv=support1(b), bad, tnv, index;
2318 ideal sub=psubset(b),M;
2319 sbv=supports(sub);
2320 nv=Nabv(hvs,av,bv);
2321 M=idMaken(mv);
2322 index = gensindex(M, idsrRing(h));
2323 n=nv.size();
2324 ring r=currRing;
2325 if(n > 0)
2326 {
2327 tnv=tnab(hvs,nv,sbv);
2328 for(i=0;i<tnv.size();i++)
2329 {
2330 co=tnv[i];
2331 bad.push_back(co+1);
2332 }
2333 for(i=0;i<n;i++)
2334 {
2335 for(j=i+1;j<n;j++)
2336 {
2337 if(nabtconditionv(hvs,nv[i],nv[j]))
2338 {
2339 good=listsinsertlist(good,i+1,j+1);
2340 }
2341 else
2342 {
2343 ;
2344 }
2345 }
2346 }
2347 solve=eli2(n,bad,good);
2348 if(bv.size()!=1)
2349 {;
2350 //PrintS("This is the solution of coefficients:\n");
2351 //listsprint(solve);
2352 }
2353 else
2354 {
2355 std::vector<int> su=make1(n);
2356 std::vector<std::vector<int> > suu;
2357 suu.push_back(su);
2358 equmab(n);
2359 solve=vecqring(solve,suu);
2360 //PrintS("This is the solution of coefficients:\n");
2361 //listsprint(solve);
2362 rChangeCurrRing(r);
2363 }
2364 solve=value1(mv,nv,solve,av,bv);
2365 }
2366 else
2367 {
2368 //PrintS("No element considered here!\n");
2369 solve.clear();
2370 }
2371 //PrintS("This is the solution of final coefficients:\n");
2372 //listsprint(solve);
2374 intvec *sl=Tmat(solve);
2375 //sl->show(0,0);
2376 return sl;
2377}
2378
2379//for debugging
2380static void T1(ideal h)
2381{
2382 ideal bi=findb(h),ai;
2383 int mm=0;
2384 id_print(bi);
2385 poly a,b;
2386 std::vector<std::vector<int> > solve;
2387 for(int i=0;i<IDELEMS(bi);i++)
2388 {
2389 //PrintS("This is aset according to:");
2390 b=pCopy(bi->m[i]);
2391 pWrite(b);
2392 ai=finda(h,b,0);
2393 if(!idIs0(ai))
2394 {
2395 id_print(ai);
2396 for(int j=0;j<IDELEMS(ai);j++)
2397 {
2398 //PrintS("This is a:");
2399 a=pCopy(ai->m[j]);
2400 //pWrite(a);
2401 intvec * solve=gradedpiece1n(h, a, b);
2402 if (IMATELEM(*solve,1,1)!=10)
2403 mm++;
2404 }
2405 }
2406 }
2407 Print("Finished %d!\n",mm);
2408}
2409
2410static bool condition2for2nv(std::vector<std::vector<int> > hvs, std::vector<int> pv, std::vector<int> qv, std::vector<int> fv)
2411{
2412 std::vector<int> vec=vecUnion(pv,qv);
2413 vec=vecUnion(vec,fv);
2414 if(vInvsl(vec,hvs))
2415 {
2416 //PrintS("condition2for2 yes\n");
2417 return (true);
2418 }
2419 //PrintS("condition2for2 no\n");
2420 return (false);
2421}
2422
2423//for subspace of T2(find all the possible values of alpha)
2424static std::vector<int> findalphan(std::vector<std::vector<int> > N, std::vector<int> tN)
2425{
2426 int i;std::vector<int> alset,vec;
2427 for(i=0;i<N.size();i++)
2428 {
2429 // vec=N[i];
2430 if(!IsinL(i,tN))
2431 {
2432 alset.push_back(i);
2433 }
2434 }
2435 //listprint(alset);
2436 return alset;
2437}
2438
2439//subspace of T^2 (nab method)
2440static std::vector<std::vector<int> > subspacetn(std::vector<std::vector<int> > N, std::vector<int> tN, std::vector<std::vector<int> > ntvs)
2441{
2442 int i;
2443 std::vector<int> alset=findalphan(N,tN), subase;
2444 std::vector<std::vector<int> > subases;
2445 for(i=0;i<alset.size();i++)
2446 {
2447 subase=subspacet1(alset[i],ntvs);
2448 subases.push_back(subase);
2449 }
2450 //PrintS("These are the bases for the subspace:\n");
2451 //listsprint(subases);
2452 return subases;
2453}
2454
2455//mts Mabt
2456//nts Nabt
2457//mvs Mab
2458//nvs Nab
2459static std::vector<std::vector<int> > value2(std::vector<std::vector<int> > mvs, std::vector<std::vector<int> > nvs, std::vector<std::vector<int> > mts, std::vector<std::vector<int> > nts, std::vector<std::vector<int> > vecs,std::vector<int> av, std::vector<int> bv)
2460{
2461 int row,col,j;
2462 std::vector<int> pv,qv, base;
2463 std::vector<std::vector<int> > bases;
2464 //PrintS("This is the nabt:\n");
2465 //listsprint(nts);
2466 //PrintS("nabt ends:\n");
2467 //PrintS("This is the mabt:\n");
2468 //listsprint(mts);
2469 //PrintS("mabt ends:\n");
2470 for(unsigned t=0;t<vecs.size();t++)
2471 {
2472 for(unsigned i=0;i<mts.size();i++)
2473 {
2474 row=mts[i][0];
2475 col=mts[i][1];
2476 pv=phimage(mvs[row],av,bv);
2477 qv=phimage(mvs[col],av,bv);
2478 if(vEvl(pv,qv))
2479 base.push_back(0);
2480 else
2481 {
2482 for(j=0;j<nts.size();j++)
2483 {
2484 row=nts[j][0];
2485 col=nts[j][1];
2486 if(vEvl(pv,nvs[row])&&vEvl(qv,nvs[col]))
2487 {
2488 base.push_back(vecs[t][j]);break;
2489 }
2490 else if(vEvl(pv,nvs[col])&&vEvl(qv,nvs[row]))
2491 {
2492 base.push_back(-vecs[t][j]);break;
2493 }
2494 }
2495 if(j==nts.size()) {base.push_back(0);}
2496 }
2497 }
2498 if(base.size()!=mts.size())
2499 {
2500 WerrorS("Errors in Values Finding(value2)!");
2501 //WerrorS("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1");
2502 usleep(1000000);
2503 assert(false);
2504 }
2505 bases.push_back(base);
2506 base.clear();
2507 }
2508 return bases;
2509}
2510
2511static ideal genst(ideal h, poly a, poly b)
2512{
2513 std::vector<std::vector<int> > hvs=supports(h),mv,mts;
2514 std::vector<int> av=support1(a), bv=support1(b);
2515 mv=Mabv(h,a,b);
2516 mts=mabtv(hvs,mv,av,bv);
2517 std::vector<std::vector<poly> > pvs=idMakei(mv,mts);
2518 ideal gens=idInit(1,1);
2519 for(unsigned i=0;i<pvs.size();i++)
2520 {
2521 idInsertPoly(gens,pvs[i][0]);
2522 idInsertPoly(gens,pvs[i][1]);
2523 }
2524 idSkipZeroes(gens);
2525 return (gens);
2526}
2527
2528static intvec * gradedpiece2n(ideal h,poly a,poly b)
2529{
2530 int i,j,t,n;
2531 std::vector<std::vector<int> > hvs=supports(h),nv,mv,mts,sbv,vecs,vars,ntvs,solve;
2532 std::vector<int> av=support1(a), bv=support1(b),tnv,vec,var;
2533 ideal sub=psubset(b);
2534 sbv=supports(sub);
2535 nv=Nabv(hvs,av,bv);
2536 n=nv.size();
2537 tnv=tnab(hvs,nv,sbv);
2538 ring r=currRing;
2539 mv=Mabv(h,a,b);
2540 mts=mabtv(hvs,mv,av,bv);
2541 //PrintS("The relations are:\n");
2542 //listsprint(mts);
2543 //PrintS("The homomorphism should map onto:\n");
2544 //lpsprint(idMakei(mv,mts));
2545 if(n>0)
2546 {
2547 ntvs=nabtv( hvs, nv, av, bv);
2548 //PrintS("The current homomorphism map onto###:\n");
2549 //lpsprint(idMakei(nv,ntvs));
2550 int l=ntvs.size();
2551 for(int t0=0;t0<l;t0++)
2552 {
2553 i=ntvs[t0][0];
2554 j=ntvs[t0][1];
2555 if(tNab(hvs,nv[i],sbv)&&tNab(hvs,nv[j],sbv))//condition 1
2556 {
2557 vec=makeequation(t0+1,0,0);
2558 vecs.push_back(vec);
2559 vec.clear();
2560 }
2561 for(int t1=t0+1;t1<ntvs.size();t1++)
2562 {
2563 for(int t2=t1+1;t2<ntvs.size();t2++)
2564 {
2565 if(ntvs[t0][0]==ntvs[t1][0]&&ntvs[t1][1]==ntvs[t2][1]&&ntvs[t0][1]==ntvs[t2][0])
2566 {
2567 i=ntvs[t0][0];
2568 j=ntvs[t0][1];
2569 t=ntvs[t1][1];
2570 if(condition2for2nv(hvs,nv[i],nv[j],nv[t]))
2571 {
2572 vec=makeequation(t0+1,t1+1,t2+1);
2573 vecs.push_back(vec);
2574 vec.clear();
2575 }
2576 }
2577 }
2578 }
2579 }
2580 //PrintS("this is EQUATIONS:\n");
2581 //listsprint(vecs);
2582 if(n==1) l=1;
2583 equmab(l);
2584 ideal id_re=idMake3(vecs);
2585 //id_print(id_re);
2586 std::vector<std::vector<int> > re=getvector(id_re,l);
2587 //PrintS("this is the solution for ideal :\n");
2588 //listsprint(re);
2589 rChangeCurrRing(r);
2590 std::vector<std::vector<int> > sub=subspacetn(nv, tnv,ntvs);
2591 //PrintS("this is the solution for subspace:\n");
2592 //listsprint(sub);
2593 equmab(l);
2594 solve=vecqring(re, sub);
2595 //PrintS("This is the solution of coefficients:\n");
2596 //listsprint(solve);
2597 rChangeCurrRing(r);
2598 solve=value2(mv,nv,mts,ntvs,solve,av,bv);
2599 }
2600 else
2601 solve.clear();
2602 intvec *sl=Tmat(solve);
2603 return sl;
2604}
2605
2606//for debugging
2607static void T2(ideal h)
2608{
2609 ideal bi=findb(h),ai;
2610 id_print(bi);
2611 poly a,b;
2612 int mm=0,gp=0;
2613 std::vector<int> bv,av;
2614 std::vector<std::vector<int> > solve;
2615 for(int i=0;i<IDELEMS(bi);i++)
2616 {
2617 b=pCopy(bi->m[i]);
2618 //bv=support1(b);
2619 //PrintS("This is aset according to:");
2620 pWrite(b);
2621//if(bv.size()==2)
2622 //{
2623 ai=finda(h,b,0);
2624 if(!idIs0(ai))
2625 {
2626 PrintS("This is a set according to current b:\n");
2627 id_print(ai);
2628 for(int j=0;j<IDELEMS(ai);j++)
2629 {
2630 PrintS("This is a:");
2631 a=pCopy(ai->m[j]);
2632 pWrite(a);
2633 PrintS("This is b:");
2634 pWrite(b);
2636 delete solve;
2637 gp++;
2638 }
2639 }
2640 mm=mm+1;
2641 }
2642 if(mm==IDELEMS(bi))
2643 PrintS("Finished!\n");
2644 Print("There are %d graded pieces in total.\n",gp);
2645}
2646
2647/*****************************for links*******************************************/
2648//the image phi(pv)=pv minus av minus bv
2649static std::vector<int> phimagel(std::vector<int> fv, std::vector<int> av, std::vector<int> bv)
2650{
2651 std::vector<int> nv;
2652 nv=vecMinus(fv,bv);
2653 nv=vecMinus(nv,av);
2654 return nv;
2655}
2656
2657//mvs and nvs are the supports of ideal Mab and Nab
2658//vecs is the solution of nab
2659static std::vector<std::vector<int> > value1l(std::vector<std::vector<int> > mvs, std::vector<std::vector<int> > lks, std::vector<std::vector<int> > vecs,std::vector<int> av, std::vector<int> bv)
2660{
2661 int j;
2662 std::vector<int> pv;
2663 std::vector<int> base;
2664 std::vector<std::vector<int> > bases;
2665 for(unsigned t=0;t<vecs.size();t++)
2666 {
2667 for(unsigned i=0;i<mvs.size();i++)
2668 {
2669 pv=phimagel(mvs[i], av, bv);
2670 for(j=0;j<lks.size();j++)
2671 {
2672 if(vEvl(pv,lks[j]))
2673 {
2674 base.push_back(vecs[t][j]);break;
2675 }
2676 }
2677 //if(j==lks.size()) {base.push_back(0);}
2678 }
2679 if(base.size()!=mvs.size())
2680 {
2681 WerrorS("Errors in Values Finding(value1l)!");
2682 usleep(1000000);
2683 assert(false);
2684 }
2685 bases.push_back(base);
2686 base.clear();
2687 }
2688 return bases;
2689}
2690
2691/***************************************************/
2693/**************************************************/
2694
2695static void TimeShow(clock_t t_construct, clock_t t_solve, clock_t t_value ,clock_t t_total)
2696{
2697 Print("The time of value matching for first order deformation: %.2f sec ;\n", ((double) t_value)/CLOCKS_PER_SEC);
2698 Print("The total time of fpiece: %.2f sec ;\n", ((double) t_total)/CLOCKS_PER_SEC);
2699 Print("The time of equations construction for fpiece: %.2f sec ;\n", ((double) t_construct)/CLOCKS_PER_SEC);
2700 Print("The total time of equations solving for fpiece: %.2f sec ;\n", ((double) t_solve)/CLOCKS_PER_SEC);
2701 PrintS("__________________________________________________________\n");
2702}
2703
2704static std::vector<std::vector<int> > gpl(ideal h,poly a,poly b)
2705{
2706 int i,j,co;
2707 std::vector<std::vector<int> > hvs=supports(h),sbv,nv,mv,good,solve;
2708 std::vector<int> av=support1(a), bv=support1(b),index,bad,tnv;
2709 ideal sub=psubset(b);
2710 sbv=supports(sub);
2711 nv=Nabv(hvs,av,bv);
2712 mv=Mabv(h,a,b);
2713 ideal M=idMaken(mv);
2714 index = gensindex(M, idsrRing(h));
2715 int n=nv.size();
2716 ring r=currRing;
2717 t_begin=clock();
2718 if(n > 0)
2719 {
2720 tnv=tnab(hvs,nv,sbv);
2721 for(i=0;i<tnv.size();i++)
2722 {
2723 co=tnv[i];
2724 bad.push_back(co+1);
2725 }
2726 for(i=0;i<n;i++)
2727 {
2728 for(j=i+1;j<n;j++)
2729 {
2730 if(nabtconditionv(hvs,nv[i],nv[j]))
2731 {
2732 good=listsinsertlist(good,i+1,j+1);
2733 }
2734 else
2735 {
2736 ;
2737 }
2738 }
2739 }
2741 t_begin=clock();
2742 solve=eli2(n,bad,good);
2743 t_solve=t_solve+clock()-t_begin;
2744 if(bv.size()!=1)
2745 {;
2746 }
2747 else
2748 {
2749 std::vector<int> su=make1(n);
2750 std::vector<std::vector<int> > suu;
2751 suu.push_back(su);
2752 equmab(n);
2753 solve=vecqring(solve,suu);
2754 rChangeCurrRing(r);
2755 }
2756 }
2757 else
2758 {
2759 solve.clear();
2760 }
2761 //listsprint(solve);
2762 //sl->show(0,0);
2763 return solve;
2764}
2765
2766//T^1
2767//only need to consider the links of a, and reduce a to empty set
2768static intvec * gradedpiece1nl(ideal h,poly a,poly b, int set)
2769{
2770 t_start=clock();
2771 poly e=pOne();
2772 std::vector<int> av=support1(a),bv=support1(b),index, em;
2773 std::vector<std::vector<int> > solve, hvs=supports(h), lks=links(a,h), mv=Mabv(h,a,b), nvl;
2774 ideal id_links=idMaken(lks);
2775 ideal M=idMaken(mv);
2776 index = gensindex(M, idsrRing(h));
2777 solve=gpl(id_links,e,b);
2778 t_mark=clock();
2779 nvl=Nabv(lks,em,bv);
2780 solve=value1l(mv, nvl , solve, av, bv);
2781 if(set==1)
2782 {
2784 }
2785 intvec *sl=Tmat(solve);
2786 t_value=t_value+clock()-t_mark;
2787 t_total=t_total+clock()-t_start;
2788 return sl;
2789}
2790
2791//for finding values of T^2
2792static std::vector<std::vector<int> > value2l(std::vector<std::vector<int> > mvs, std::vector<std::vector<int> > lks, std::vector<std::vector<int> > mts, std::vector<std::vector<int> > lkts, std::vector<std::vector<int> > vecs,std::vector<int> av, std::vector<int> bv)
2793{
2794 std::vector<int> pv,qv,base;
2795 int row,col,j;
2796 std::vector<std::vector<int> > bases;
2797 if(vecs.size()==0)
2798 {
2799
2800 }
2801 for(unsigned t=0;t<vecs.size();t++)
2802 {
2803 for(unsigned i=0;i<mts.size();i++)
2804 {
2805 row=mts[i][0];
2806 col=mts[i][1];
2807 pv=phimagel(mvs[row],av,bv);
2808 qv=phimagel(mvs[col],av,bv);
2809 if(vEvl(pv,qv))
2810 base.push_back(0);
2811 else
2812 {
2813 for(j=0;j<lkts.size();j++)
2814 {
2815 row=lkts[j][0];
2816 col=lkts[j][1];
2817 if(vEvl(pv,lks[row])&&vEvl(qv,lks[col]))
2818 {
2819 base.push_back(vecs[t][j]);break;
2820 }
2821 else if(vEvl(qv,lks[row])&&vEvl(pv,lks[col]))
2822 {
2823 base.push_back(-vecs[t][j]);break;
2824 }
2825 }
2826 //if(j==lkts.size())
2827 //{
2828 //base.push_back(0);
2829 //}
2830 }
2831 }
2832 if(base.size()!=mts.size())
2833 {
2834 WerrorS("Errors in Values Finding!");
2835 //WerrorS("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1");
2836 usleep(1000000);
2837 assert(false);
2838 }
2839 bases.push_back(base);
2840 base.clear();
2841 }
2842 return bases;
2843}
2844
2845static std::vector<std::vector<int> > gpl2(ideal h,poly a,poly b)
2846{
2847 int i,j,t,n;
2848 std::vector<std::vector<int> > hvs=supports(h),sbv,nv,mv,mts,vecs,vars,ntvs,solve;
2849 std::vector<int> av=support1(a), bv=support1(b),vec,var,tnv;
2850 ideal sub=psubset(b);
2851 sbv=supports(sub);
2852 nv=Nabv(hvs,av,bv);
2853 n=nv.size();
2854 tnv=tnab(hvs,nv,sbv);
2855 ring r=currRing;
2856 mv=Mabv(h,a,b);
2857 mts=mabtv(hvs,mv,av,bv);
2858 if(n>0)
2859 {
2860 ntvs=nabtv( hvs, nv, av, bv);
2861 int l=ntvs.size();
2862 if(l>0)
2863 {
2864 for(int t0=0;t0<l;t0++)
2865 {
2866 i=ntvs[t0][0];
2867 j=ntvs[t0][1];
2868 if(tNab(hvs,nv[i],sbv)&&tNab(hvs,nv[j],sbv))//condition 1
2869 {
2870 vec=makeequation(t0+1,0,0);
2871 vecs.push_back(vec);
2872 vec.clear();
2873 }
2874 for(int t1=t0+1;t1<ntvs.size();t1++)
2875 {
2876 for(int t2=t1+1;t2<ntvs.size();t2++)
2877 {
2878 if(ntvs[t0][0]==ntvs[t1][0]&&ntvs[t1][1]==ntvs[t2][1]&&ntvs[t0][1]==ntvs[t2][0])
2879 {
2880 i=ntvs[t0][0];
2881 j=ntvs[t0][1];
2882 t=ntvs[t1][1];
2883 if(condition2for2nv(hvs,nv[i],nv[j],nv[t]))
2884 {
2885 vec=makeequation(t0+1,t1+1,t2+1);
2886 vecs.push_back(vec);
2887 vec.clear();
2888 }
2889 }
2890 }
2891 }
2892 }
2893 if(n==1) {l=1;}
2894 equmab(l);
2895 ideal id_re=idMake3(vecs);
2896 std::vector<std::vector<int> > re=getvector(id_re,l);
2897 rChangeCurrRing(r);
2898 std::vector<std::vector<int> > sub=subspacetn(nv, tnv,ntvs);
2899 equmab(l);
2900 solve=vecqring(re, sub);
2901 rChangeCurrRing(r);
2902 }
2903 else
2904 {
2905 solve.clear();
2906 }
2907 }
2908 else
2909 solve.clear();
2910 return solve;
2911}
2912
2913static intvec * gradedpiece2nl(ideal h,poly a,poly b)
2914{
2915 poly e=pOne();
2916 std::vector<int> av=support1(a), bv=support1(b), em;
2917 std::vector<std::vector<int> > hvs=supports(h), mv=Mabv(h,a,b),mts,solve,lks,nvl,ntsl;
2918 mts=mabtv(hvs,mv,av,bv);
2919 lks=links(a,h);
2920 ideal id_links=idMaken(lks);
2921//PrintS("This is the links of a:\n"); id_print(id_links);
2922 nvl=Nabv(lks,em,bv);
2923//PrintS("This is the N set:\n"); id_print(idMaken(nvl));
2924 ntsl=nabtv(lks,nvl,em,bv);
2925//PrintS("This is N^2:\n"); listsprint(ntsl);
2926 solve=gpl2(id_links,e,b);
2927//PrintS("This is pre solution of N:\n"); listsprint(solve);
2928 if(solve.size() > 0)
2929 {
2930 solve=value2l(mv, nvl, mts, ntsl, solve, av, bv);
2931 }
2932//PrintS("This is solution of N:\n"); listsprint(solve);
2933 intvec *sl=Tmat(solve);
2934 return sl;
2935}
2936
2937//for debugging
2938/*
2939void Tlink(ideal h,poly a,poly b,int n)
2940{
2941 std::vector<std::vector<int> > hvs=supports(h);
2942 std::vector<int> av=support1(a);
2943 std::vector<int> bv=support1(b);
2944 std::vector<std::vector<int> > vec=links(a, h);
2945 PrintS("This is the links of a:\n");
2946 listsprint(vec);
2947 ideal li=idMaken(vec);
2948 PrintS("This is the links of a(ideal version):\n");
2949 id_print(li);
2950 poly p=pOne();
2951 PrintS("1************************************************\n");
2952 PrintS("This is T_1 (m):\n");
2953 gradedpiece1(li,p,b);
2954 PrintS("2************************************************\n");
2955 PrintS("This is T_2 (m):\n");
2956 gradedpiece2(li,p,b);
2957 PrintS("3************************************************\n");
2958 PrintS("This is T_1 (n):\n");
2959 gradedpiece1n(li,p,b);
2960 PrintS("4************************************************\n");
2961 PrintS("This is T_2 (n):\n");
2962 gradedpiece2n(li,p,b);
2963}
2964*/
2965
2966/******************************for triangulation***********************************/
2967
2968//returns all the faces which are triangles
2969static ideal trisets(ideal h)
2970{
2971 int i;
2972 ideal ids=idInit(1,1);
2973 std::vector<int> pv;
2974 for(i=0;i<IDELEMS(h);i++)
2975 {
2976 pv= support1(h->m[i]);
2977 if(pv.size()==3)
2978 idInsertPoly(ids, pCopy(h->m[i]));
2979 }
2980 idSkipZeroes(ids);
2981 return ids;
2982}
2983
2984// case 1 new faces
2985static std::vector<std::vector<int> > triface(poly p, int vert)
2986{
2987 std::vector<int> vec, fv=support1(p);
2988 std::vector<std::vector<int> > fvs0, fvs;
2989 vec.push_back(vert);
2990 fvs.push_back(vec);
2991 fvs0=b_subsets(fv);
2992 fvs0=vsMinusv(fvs0,fv);
2993 for(unsigned i=0;i<fvs0.size();i++)
2994 {
2995 vec=fvs0[i];
2996 vec.push_back(vert);
2997 fvs.push_back(vec);
2998 }
2999 return (fvs);
3000}
3001
3002// the size of p's support must be 3
3003//returns the new complex which is a triangulation based on the face p
3004static ideal triangulations1(ideal h, poly p, int vert)
3005{
3006 std::vector<int> vec, pv=support1(p);
3007 std::vector<std::vector<int> > vecs=supports(h),vs,vs0;
3008 vs0=triface(p,vert);
3009 vecs=vsMinusv(vecs, pv);
3010 vecs=vsUnion(vecs,vs0);
3011 //PrintS("This is the new simplicial complex according to the face \n"); pWrite(p);
3012 //PrintS("is:\n");
3013 //listsprint(vecs);
3014 ideal re=idMaken(vecs);
3015 return re;
3016}
3017
3018/*
3019static ideal triangulations1(ideal h)
3020{
3021 int i,vert=currRing->N+1;
3022 std::vector<int> vec;
3023 std::vector<std::vector<int> > vecs=supports(h),vs,vs0;
3024 for (i=0;i<vecs.size();i++)
3025 {
3026 if((vecs[i]).size()==3)
3027 {
3028 vs0=triface(vecs[i],vert);
3029 vs=vsMinusv(vecs,vecs[i]);
3030 vs=vsUnion(vs,vs0);
3031 PrintS("This is the new simplicial complex according to the face \n");listprint(vecs[i]);
3032 PrintS("is:\n");
3033 listsprint(vs);
3034 }
3035 //else if((vecs[i]).size()==4)
3036 //tetraface(vecs[i]);
3037 }
3038 //ideal hh=idMaken(vs);
3039 return h;
3040}*/
3041
3042static std::vector<int> commonedge(poly p, poly q)
3043{
3044 std::vector<int> ev, fv1= support1(p), fv2= support2(q);
3045 for(unsigned i=0;i<fv1.size();i++)
3046 {
3047 if(IsinL(fv1[i], fv2))
3048 ev.push_back(fv1[i]);
3049 }
3050 return ev;
3051}
3052
3053static intvec *edgemat(poly p, poly q)
3054{
3055 intvec *m;
3056 int i;
3057 std::vector<int> dg=commonedge(p, q);
3058 int lg=dg.size();
3059 m=new intvec(lg);
3060 if(lg!=0)
3061 {
3062 m=new intvec(lg);
3063 for(i=0;i<lg;i++)
3064 {
3065 (*m)[i]=dg[i];
3066 }
3067 }
3068 return (m);
3069}
3070
3071// case 2 the new face
3072static std::vector<std::vector<int> > tetraface(poly p, poly q, int vert)
3073{
3074 std::vector<int> ev=commonedge(p, q), vec, fv1=support1(p), fv2=support1(q);
3075 std::vector<std::vector<int> > fvs1, fvs2, fvs;
3076 vec.push_back(vert);
3077 fvs.push_back(vec);
3078 fvs1=b_subsets(fv1);
3079 fvs2=b_subsets(fv2);
3080 fvs1=vsMinusv(fvs1, fv1);
3081 fvs2=vsMinusv(fvs2, fv2);
3082 fvs2=vsUnion(fvs1, fvs2);
3083 fvs2=vsMinusv(fvs2, ev);
3084 for(unsigned i=0;i<fvs2.size();i++)
3085 {
3086 vec=fvs2[i];
3087 vec.push_back(vert);
3088 fvs.push_back(vec);
3089 }
3090 return (fvs);
3091}
3092
3093//if p and q have a common edge
3094static ideal triangulations2(ideal h, poly p, poly q, int vert)
3095{
3096 std::vector<int> ev, fv1=support1(p), fv2=support1(q);
3097 std::vector<std::vector<int> > vecs=supports(h), vs1;
3098 ev=commonedge(p, q);
3099 vecs=vsMinusv(vecs, ev);
3100 vecs=vsMinusv(vecs,fv1);
3101 vecs=vsMinusv(vecs,fv2);
3102 vs1=tetraface(p, q, vert);
3103 vecs=vsUnion(vecs,vs1);
3104 ideal hh=idMaken(vecs);
3105 return hh;
3106}
3107
3108// case 2 the new face
3109static std::vector<std::vector<int> > penface(poly p, poly q, poly g, int vert)
3110{
3111 int en=0;
3112 std::vector<int> ev1=commonedge(p, q), ev2=commonedge(p, g), ev3=commonedge(q, g), ind, vec, fv1=support1(p), fv2=support1(q), fv3=support1(g);
3113 std::vector<std::vector<int> > fvs1, fvs2, fvs3, fvs, evec;
3114 evec.push_back(ev1);
3115 evec.push_back(ev2);
3116 evec.push_back(ev3);
3117 for(unsigned i=0;i<evec.size();i++)
3118 {
3119 if(evec[i].size()==2)
3120 {
3121 en++;
3122 }
3123 }
3124 if(en==2)
3125 {
3126 vec.push_back(vert);
3127 fvs.push_back(vec);
3128 fvs1=b_subsets(fv1);
3129 fvs2=b_subsets(fv2);
3130 fvs3=b_subsets(fv3);
3131 fvs1=vsMinusv(fvs1, fv1);
3132 fvs2=vsMinusv(fvs2, fv2);
3133 fvs3=vsMinusv(fvs3, fv3);
3134 fvs3=vsUnion(fvs3, fvs2);
3135 fvs3=vsUnion(fvs3, fvs1);
3136 for(unsigned i=0;i<evec.size();i++)
3137 {
3138 if(evec[i].size()==2)
3139 {
3140 fvs3=vsMinusv(fvs3, evec[i]);
3141 }
3142 }
3143 for(unsigned i=0;i<fvs3.size();i++)
3144 {
3145 vec=fvs3[i];
3146 vec.push_back(vert);
3147 fvs.push_back(vec);
3148 }
3149 }
3150 return (fvs);
3151}
3152
3153static ideal triangulations3(ideal h, poly p, poly q, poly g, int vert)
3154{
3155 std::vector<int> ev1=commonedge(p, q), ev2=commonedge(p, g), ev3=commonedge(q, g), fv1=support1(p), fv2=support1(q), fv3=support1(g);
3156 std::vector<std::vector<int> > vecs=supports(h), vs1, evec;
3157 evec.push_back(ev1);
3158 evec.push_back(ev2);
3159 evec.push_back(ev3);
3160 for(unsigned i=0;i<evec.size();i++)
3161 {
3162 if(evec[i].size()==2)
3163 {
3164 vecs=vsMinusv(vecs, evec[i]);
3165 }
3166 }
3167 vecs=vsMinusv(vecs,fv1);
3168 vecs=vsMinusv(vecs,fv2);
3169 vecs=vsMinusv(vecs,fv3);
3170 vs1=penface(p, q, g, vert);
3171 vecs=vsUnion(vecs,vs1);
3172 ideal hh=idMaken(vecs);
3173 return hh;
3174}
3175
3176//returns p's valency in h
3177//p must be a vertex
3178static int valency(ideal h, poly p)
3179{
3180 int val=0;
3181 std::vector<int> ev=support1(pCopy(p));
3182 int ver=ev[0];
3183//PrintS("the vertex is :\n"); listprint(p);
3184 std::vector<std::vector<int> > vecs=supports(idCopy(h));
3185 for(unsigned i=0;i<vecs.size();i++)
3186 {
3187 if(vecs[i].size()==2 && IsinL(ver, vecs[i]))
3188 val++;
3189 }
3190 return (val);
3191}
3192
3193/*ideal triangulations2(ideal h)
3194{
3195 int i,j,vert=currRing->N+1;
3196 std::vector<int> ev;
3197 std::vector<std::vector<int> > vecs=supports(h),vs,vs0,vs1;
3198 vs0=tetrasets(h);
3199 for (i=0;i<vs0.size();i++)
3200 {
3201 for(j=i;j<vs0.size();j++)
3202 {
3203 ev=commonedge(vs0[i],vs0[j]);
3204 if(ev.size()==2)
3205 {
3206 vecs=vsMinusv(vecs, ev);
3207 vs=vsMinusv(vecs,vs0[i]);
3208 vs=vsMinusv(vecs,vs0[j]);
3209 vs1=tetraface(vs0[i],vs0[j],vert);
3210 vs=vsUnion(vs,vs1);
3211 PrintS("This is the new simplicial complex according to the face 1 \n");listprint(vecs[i]);
3212PrintS("face 2: \n");
3213 PrintS("is:\n");
3214 listsprint(vs);
3215 }
3216 }
3217
3218 //else if((vecs[i]).size()==4)
3219 //tetraface(vecs[i]);
3220 }
3221 //ideal hh=idMaken(vs);
3222 return h;
3223}*/
3224
3225/*********************************For computation of X_n***********************************/
3226static std::vector<std::vector<int> > vsMinusvs(std::vector<std::vector<int> > vs1, std::vector<std::vector<int> > vs2)
3227{
3228 std::vector<std::vector<int> > vs=vs1;
3229 for(unsigned i=0;i<vs2.size();i++)
3230 {
3231 vs=vsMinusv(vs, vs2[i]);
3232 }
3233 return vs;
3234}
3235
3236static std::vector<std::vector<int> > vs_subsets(std::vector<std::vector<int> > vs)
3237{
3238 std::vector<std::vector<int> > sset, bv;
3239 for(unsigned i=0;i<vs.size();i++)
3240 {
3241 bv=b_subsets(vs[i]);
3242 sset=vsUnion(sset, bv);
3243 }
3244 return sset;
3245}
3246
3247static std::vector<std::vector<int> > p_constant(ideal Xo, ideal Sigma)
3248{
3249 std::vector<std::vector<int> > xs=supports(idCopy(Xo)), ss=supports(idCopy(Sigma)), fvs1;
3250 fvs1=vs_subsets(ss);
3251 fvs1=vsMinusvs(xs, fvs1);
3252 return fvs1;
3253}
3254
3255static std::vector<std::vector<int> > p_change(ideal Sigma)
3256{
3257 std::vector<std::vector<int> > ss=supports(idCopy(Sigma)), fvs;
3258 fvs=vs_subsets(ss);
3259 return (fvs);
3260}
3261
3262static std::vector<std::vector<int> > p_new(ideal Xo, ideal Sigma)
3263{
3264 int vert=0;
3265 std::vector<std::vector<int> > ss=supports(idCopy(Sigma)), fvs;
3266 for(int i=1;i<=currRing->N;i++)
3267 {
3268 for(int j=0;j<IDELEMS(Xo);j++)
3269 {
3270 if(pGetExp(Xo->m[j],i)>0)
3271 {
3272 vert=i+1;
3273 break;
3274 }
3275 }
3276 }
3277 int typ=ss.size();
3278 if(typ==1)
3279 {
3280 fvs=triface(Sigma->m[0], vert);
3281 }
3282 else if(typ==2)
3283 {
3284 fvs=tetraface(Sigma->m[0], Sigma->m[1], vert);
3285 }
3286 else
3287 {
3288 fvs=penface(Sigma->m[0], Sigma->m[1], Sigma->m[2], vert);
3289 }
3290 return (fvs);
3291}
3292
3293static ideal c_New(ideal Io, ideal sig)
3294{
3295 std::vector<std::vector<int> > vs1=p_constant(Io, sig), vs2=p_change(sig), vs3=p_new(Io, sig), vsig=supports(sig), vs;
3296 std::vector<int> ev;
3297 int ednum=vsig.size();
3298 if(ednum==2)
3299 {
3300 vsig.push_back(commonedge(sig->m[0], sig->m[1]));
3301 }
3302 else if(ednum==3)
3303 {
3304 for(int i=0;i<IDELEMS(sig);i++)
3305 {
3306 for(int j=i+1;j<IDELEMS(sig);j++)
3307 {
3308 ev=commonedge(sig->m[i], sig->m[j]);
3309 if(ev.size()==2)
3310 {
3311 vsig.push_back(ev);
3312 }
3313 }
3314 }
3315 }
3316//PrintS("the first part is:\n");id_print(idMaken(vs1));
3317//PrintS("the second part is:\n");id_print(idMaken(vsig));
3318//PrintS("the third part is:\n");id_print(idMaken(vs3));
3319 vs2=vsMinusvs(vs2, vsig);
3320//PrintS("the constant part2 is:\n");id_print(idMaken(vs2));
3321 vs=vsUnion(vs2, vs1);
3322//PrintS("the constant part is:\n");id_print(idMaken(vs));
3323 vs=vsUnion(vs, vs3);
3324//PrintS("the whole part is:\n");id_print(idMaken(vs));
3325 return(idMaken(vs));
3326}
3327
3328static std::vector<std::vector<int> > phi1(poly a, ideal Sigma)
3329{
3330 std::vector<std::vector<int> > ss=supports(idCopy(Sigma)), fvs;
3331 std::vector<int> av=support1(a), intvec, vv;
3332 for(unsigned i=0;i<ss.size();i++)
3333 {
3334 intvec=vecIntersection(ss[i], av);
3335 if(intvec.size()==av.size())
3336 {
3337 vv=vecMinus(ss[i], av);
3338 fvs.push_back(vv);
3339 }
3340 }
3341 return fvs;
3342}
3343
3344static std::vector<std::vector<int> > phi2(poly a, ideal Xo, ideal Sigma)
3345{
3346
3347 std::vector<std::vector<int> > ss=p_new(Sigma, Xo), fvs;
3348 std::vector<int> av=support1(a), intvec, vv;
3349 for(unsigned i=0;i<ss.size();i++)
3350 {
3351 intvec=vecIntersection(ss[i], av);
3352 if(intvec.size()==av.size())
3353 {
3354 vv=vecMinus(ss[i], av);
3355 fvs.push_back(vv);
3356 }
3357 }
3358 return fvs;
3359}
3360
3361static std::vector<std::vector<int> > links_new(poly a, ideal Xo, ideal Sigma, int vert, int ord)
3362{
3363 std::vector<int> av=support1(a);
3364 std::vector<std::vector<int> > lko, lkn, lk1, lk2;
3365 lko=links(a, Xo);
3366 if(ord==1)
3367 return lko;
3368 if(ord==2)
3369 {
3370 lk1=phi1(a, Sigma);
3371 lk2=phi2(a, Xo, Sigma);
3372 lkn=vsMinusvs(lko, lk1);
3373 lkn=vsUnion(lkn, lk2);
3374 return lkn;
3375 }
3376 if(ord==3)
3377 {
3378 lkn=phi2(a, Xo, Sigma);
3379 return lkn;
3380 }
3381 WerrorS("Cannot find the links smartly!");
3382 return lko;
3383}
3384
3385//returns 1 if there is a real divisor of b not in Xs
3386static int existIn(poly b, ideal Xs)
3387{
3388 std::vector<int> bv=support1(pCopy(b));
3389 std::vector<std::vector<int> > xvs=supports(idCopy(Xs)), bs=b_subsets(bv);
3390 bs=vsMinusv(bs, bv);
3391 for(unsigned i=0;i<bs.size();i++)
3392 {
3393 if(!vInvsl(bs[i], xvs))
3394 {
3395 return 1;
3396 }
3397 }
3398 return 0;
3399}
3400
3401static int isoNum(poly p, ideal I, poly a, poly b)
3402{
3403 int i;
3404 std::vector<std::vector<int> > vs=supports(idCopy(I));
3405 std::vector<int> v1=support1(a), v2=support1(b), v=support1(p);
3406 std::vector<int> vp, iv=phimagel(v, v1, v2);
3407 for(i=0;i<IDELEMS(I);i++)
3408 {
3409 vp=support1(pCopy(I->m[i]));
3410 if(vEvl(iv, phimagel(vp, v1, v2)))
3411 {
3412 return (i+1);
3413 }
3414 }
3415 return (0);
3416}
3417
3418static int ifIso(poly p, poly q, poly f, poly g, poly a, poly b)
3419{
3420 std::vector<int> va=support1(a), vb=support1(b), vp=support1(p), vq=support1(q), vf=support1(f), vg=support1(g);
3421 std::vector<int> v1=phimagel(vp, va, vb), v2=phimagel(vq, va, vb), v3=phimagel(vf, va, vb), v4=phimagel(vg, va, vb);
3422 if((vEvl(v1, v3)&& vEvl(v2,v4))||(vEvl(v1, v4)&& vEvl(v2,v3)) )
3423 {
3424 return (1);
3425 }
3426 return (0);
3427}
3428
3429static ideal idMinusp(ideal I, poly p)
3430{
3431 ideal h=idInit(1,1);
3432 int i;
3433 for(i=0;i<IDELEMS(I);i++)
3434 {
3435 if(!p_EqualPolys(I->m[i], p, currRing))
3436 {
3437 idInsertPoly(h, pCopy(I->m[i]));
3438 }
3439 }
3440 idSkipZeroes(h);
3441 return h;
3442}
3443
3444/****************************for the interface of .lib*********************************/
3445
3446static std::vector<int> v_minus(std::vector<int> v1, std::vector<int> v2)
3447{
3448 std::vector<int> vec;
3449 for(unsigned i=0;i<v1.size();i++)
3450 {
3451 vec.push_back(v1[i]-v2[i]);
3452 }
3453 return vec;
3454}
3455
3456static std::vector<int> gdegree(poly a, poly b)
3457{
3458 int i;
3459 std::vector<int> av,bv;
3460 for(i=1;i<=currRing->N;i++)
3461 {
3462 av.push_back(pGetExp(a,i));
3463 bv.push_back(pGetExp(b,i));
3464 }
3465 std::vector<int> vec=v_minus(av,bv);
3466 //PrintS("The degree is:\n");
3467 //listprint(vec);
3468 return vec;
3469}
3470
3471/********************************for stellar subdivision******************************/
3472
3473static std::vector<std::vector<int> > star(poly a, ideal h)
3474{
3475 int i;
3476 std::vector<std::vector<int> > st,X=supports(h);
3477 std::vector<int> U,av=support1(a);
3478 for(i=0;i<X.size();i++)
3479 {
3480 U=vecUnion(av,X[i]);
3481 if(vInvsl(U,X))
3482 {
3483 st.push_back(X[i]);
3484 }
3485 }
3486 return st;
3487}
3488
3489static std::vector<std::vector<int> > boundary(poly a)
3490{
3491 std::vector<int> av=support1(a), vec;
3492 std::vector<std::vector<int> > vecs;
3493 vecs=b_subsets(av);
3494 vecs.push_back(vec);
3495 vecs=vsMinusv(vecs, av);
3496 return vecs;
3497}
3498
3499static std::vector<std::vector<int> > stellarsub(poly a, ideal h)
3500{
3501 std::vector<std::vector<int> > vecs_minus, vecs_plus, lk=links(a,h), hvs=supports(h), sub, bys=boundary(a);
3502 std::vector<int> av=support1(a), vec, vec_n;
3503 int i,j,vert=0;
3504 for(i=1;i<=currRing->N;i++)
3505 {
3506 for(j=0;j<IDELEMS(h);j++)
3507 {
3508 if(pGetExp(h->m[j],i)>0)
3509 {
3510 vert=i+1;
3511 break;
3512 }
3513 }
3514 }
3515 vec_n.push_back(vert);
3516 for(i=0;i<lk.size();i++)
3517 {
3518 vec=vecUnion(av, lk[i]);
3519 vecs_minus.push_back(vec);
3520 for(j=0;j<bys.size();j++)
3521 {
3522 vec=vecUnion(lk[i], vec_n);
3523 vec=vecUnion(vec, bys[j]);
3524 vecs_plus.push_back(vec);
3525 }
3526 }
3527 sub=vsMinusvs(hvs, vecs_minus);
3528 sub=vsUnion(sub, vecs_plus);
3529 return(sub);
3530}
3531
3532static std::vector<std::vector<int> > bsubsets_1(poly b)
3533{
3534 std::vector<int> bvs=support1(b), vs;
3535 std::vector<std::vector<int> > bset;
3536 for(unsigned i=0;i<bvs.size();i++)
3537 {
3538 for(int j=0;j!=i; j++)
3539 {
3540 vs.push_back(bvs[j]);
3541 }
3542 bset.push_back(vs);
3543 vs.resize(0);
3544 }
3545 return bset;
3546}
3547
3548/***************************for time testing******************************/
3549static ideal T_1h(ideal h)
3550{
3551 int i, j;
3552 //std::vector < intvec > T1;
3553 ideal ai=p_a(h), bi;
3554 //intvec *L;
3555 for(i=0;i<IDELEMS(ai);i++)
3556 {
3557 bi=p_b(h,ai->m[i]);
3558 if(!idIs0(bi))
3559 {
3560 for(j=0;j<IDELEMS(bi);j++)
3561 {
3562 //PrintS("This is for:\n");pWrite(ai->m[i]); pWrite(bi->m[j]);
3563 gradedpiece1nl(h,ai->m[i],bi->m[j], 0);
3564 //PrintS("Succeed!\n");
3565 //T1.push_back(L);
3566 }
3567 }
3568 }
3570 return h;
3571}
3572
3573/**************************************interface T1****************************************/
3574/*
3575static BOOLEAN makeqring(leftv res, leftv args)
3576{
3577 leftv h=args;
3578 ideal h2= id_complement( hh);
3579 if((h != NULL)&&(h->Typ() == POLY_CMD))
3580 {
3581 poly p= (poly)h->Data();
3582 h = h->next;
3583 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3584 {
3585 ideal hh=(ideal)h->Data();
3586 ideal h2=id_complement(hh);
3587 ideal h1=id_Init(1,1);
3588 idInsertPoly(h1,p);
3589 ideal gb=kStd(h2,NULL,testHomog,NULL,NULL,0,0,NULL);
3590 ideal idq=kNF(gb,NULL,h1);
3591 idSkipZeroes(h1);
3592 res->rtyp =POLY_CMD;
3593 res->data =h1->m[0];
3594 }
3595 }
3596 }
3597 return false;
3598}*/
3599
3601{
3602 leftv h=args;
3603 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3604 {
3605 ideal hh=(ideal)h->Data();
3606 res->rtyp =IDEAL_CMD;
3607 res->data =idsrRing(hh);
3608 }
3609 return false;
3610}
3611
3613{
3614 leftv h=args;
3615 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3616 {
3617 ideal hh=(ideal)h->Data();
3618 ideal h2= id_complement(hh);
3619 res->rtyp =IDEAL_CMD;
3620 res->data =h2;
3621 }
3622 return false;
3623}
3624
3626{
3627 leftv h=args;
3628 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3629 {
3630 ideal hh=(ideal)h->Data();
3631 res->rtyp =IDEAL_CMD;
3632 res->data =T_1h(hh);
3633 }
3634 return false;
3635}
3636
3638{
3639 leftv h=args;
3640 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3641 {
3642 ideal h1= (ideal)h->Data();
3643 h = h->next;
3644 if((h != NULL)&&(h->Typ() == POLY_CMD))
3645 {
3646 poly p= (poly)h->Data();
3647 h = h->next;
3648 if((h != NULL)&&(h->Typ() == POLY_CMD))
3649 {
3650 poly q= (poly)h->Data();
3651 res->rtyp =IDEAL_CMD;
3652 res->data =mingens(h1,p,q);
3653 }
3654 }
3655 }
3656 return false;
3657}
3658
3659static intvec *dmat(poly a, poly b)
3660{
3661 intvec *m;
3662 int i;
3663 std::vector<int> dg=gdegree(a,b);
3664 int lg=dg.size();
3665 m=new intvec(lg);
3666 if(lg!=0)
3667 {
3668 m=new intvec(lg);
3669 for(i=0;i<lg;i++)
3670 {
3671 (*m)[i]=dg[i];
3672 }
3673 }
3674 return (m);
3675}
3676
3677static BOOLEAN gd(leftv res, leftv args)
3678{
3679 leftv h=args;
3680 if((h != NULL)&&(h->Typ() == POLY_CMD))
3681 {
3682 poly p= (poly)h->Data();
3683 h = h->next;
3684 if((h != NULL)&&(h->Typ() == POLY_CMD))
3685 {
3686 poly q= (poly)h->Data();
3687 res->rtyp =INTVEC_CMD;
3688 res->data =dmat(p,q);
3689 }
3690 }
3691 return false;
3692}
3693
3695{
3696 leftv h=args;
3697 if((h != NULL)&&(h->Typ() == POLY_CMD))
3698 {
3699 poly p= (poly)h->Data();
3700 h = h->next;
3701 if((h != NULL)&&(h->Typ() == POLY_CMD))
3702 {
3703 poly q= (poly)h->Data();
3704 res->rtyp =INTVEC_CMD;
3705 res->data =edgemat(p,q);
3706 }
3707 }
3708 return false;
3709}
3710
3711static BOOLEAN fb(leftv res, leftv args)
3712{
3713 leftv h=args;
3714 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3715 {
3716 ideal h1= (ideal)h->Data();
3717 res->rtyp =IDEAL_CMD;
3718 res->data =findb(h1);
3719 }
3720 return false;
3721}
3722
3723static BOOLEAN pa(leftv res, leftv args)
3724{
3725 leftv h=args;
3726 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3727 {
3728 ideal h1= (ideal)h->Data();
3729 res->rtyp =IDEAL_CMD;
3730 res->data =p_a(h1);
3731 }
3732 return false;
3733}
3734
3736{
3737 leftv h=args;
3738 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3739 {
3740 ideal h1= (ideal)h->Data();
3741 res->rtyp =IDEAL_CMD;
3742 res->data =complementsimplex(h1);
3743 }
3744 return false;
3745}
3746
3747static BOOLEAN pb(leftv res, leftv args)
3748{
3749 leftv h=args;
3750 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3751 {
3752 ideal h1= (ideal)h->Data();
3753 h = h->next;
3754 if((h != NULL)&&(h->Typ() == POLY_CMD))
3755 {
3756 poly p= (poly)h->Data();
3757 res->rtyp =IDEAL_CMD;
3758 res->data =p_b(h1,p);
3759 }
3760 }
3761 return false;
3762}
3763
3764static BOOLEAN fa(leftv res, leftv args)
3765{
3766 leftv h=args;
3767 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3768 {
3769 ideal h1= (ideal)h->Data();
3770 h = h->next;
3771 if((h != NULL)&&(h->Typ() == POLY_CMD))
3772 {
3773 poly q= (poly)h->Data();
3774 h = h->next;
3775 if((h != NULL)&&(h->Typ() == INT_CMD))
3776 {
3777 int d= (int)(long)h->Data();
3778 res->rtyp =IDEAL_CMD;
3779 res->data =finda(h1,q,d);
3780 }
3781 }
3782 }
3783 return false;
3784}
3785
3787{
3788 leftv h=args;
3789 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3790 {
3791 ideal h1= (ideal)h->Data();
3792 h = h->next;
3793 if((h != NULL)&&(h->Typ() == POLY_CMD))
3794 {
3795 poly p= (poly)h->Data();
3796 h = h->next;
3797 if((h != NULL)&&(h->Typ() == POLY_CMD))
3798 {
3799 poly q= (poly)h->Data();
3800 res->rtyp =INTVEC_CMD;
3801 res->data =gradedpiece1n(h1,p,q);
3802 }
3803 }
3804 }
3805 return false;
3806}
3807
3809{
3810 leftv h=args;
3811 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3812 {
3813 ideal h1= (ideal)h->Data();
3814 h = h->next;
3815 if((h != NULL)&&(h->Typ() == POLY_CMD))
3816 {
3817 poly p= (poly)h->Data();
3818 h = h->next;
3819 if((h != NULL)&&(h->Typ() == POLY_CMD))
3820 {
3821 poly q= (poly)h->Data();
3822 h = h->next;
3823 if((h != NULL)&&(h->Typ() == INT_CMD))
3824 {
3825 int d= (int)(long)h->Data();
3826 res->rtyp =INTVEC_CMD;
3827 res->data =gradedpiece1nl(h1,p,q,d);
3828 }
3829 }
3830 }
3831 }
3832 return false;
3833}
3834
3836{
3837 leftv h=args;
3838 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3839 {
3840 ideal h1= (ideal)h->Data();
3841 h = h->next;
3842 if((h != NULL)&&(h->Typ() == POLY_CMD))
3843 {
3844 poly p= (poly)h->Data();
3845 h = h->next;
3846 if((h != NULL)&&(h->Typ() == POLY_CMD))
3847 {
3848 poly q= (poly)h->Data();
3849 res->rtyp =IDEAL_CMD;
3850 res->data =genst(h1,p,q);
3851 }
3852 }
3853 }
3854 return false;
3855}
3856
3858{
3859 leftv h=args;
3860 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3861 {
3862 ideal h1= (ideal)h->Data();
3863 h = h->next;
3864 if((h != NULL)&&(h->Typ() == POLY_CMD))
3865 {
3866 poly p= (poly)h->Data();
3867 h = h->next;
3868 if((h != NULL)&&(h->Typ() == POLY_CMD))
3869 {
3870 poly q= (poly)h->Data();
3871 res->rtyp =INTVEC_CMD;
3872 res->data =gradedpiece2n(h1,p,q);
3873 }
3874 }
3875 }
3876 return false;
3877}
3878
3880{
3881 leftv h=args;
3882 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3883 {
3884 ideal h1= (ideal)h->Data();
3885 h = h->next;
3886 if((h != NULL)&&(h->Typ() == POLY_CMD))
3887 {
3888 poly p= (poly)h->Data();
3889 h = h->next;
3890 if((h != NULL)&&(h->Typ() == POLY_CMD))
3891 {
3892 poly q= (poly)h->Data();
3893 res->rtyp =INTVEC_CMD;
3894 res->data =gradedpiece2nl(h1,p,q);
3895 }
3896 }
3897 }
3898 return false;
3899}
3900
3902{
3903 leftv h=args;
3904 if((h != NULL)&&(h->Typ() == POLY_CMD))
3905 {
3906 poly p= (poly)h->Data();
3907 h = h->next;
3908 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3909 {
3910 ideal h1= (ideal)h->Data();
3911 res->rtyp =IDEAL_CMD;
3912 std::vector<std::vector<int> > vecs=links(p,h1);
3913 res->data =idMaken(vecs);
3914 }
3915 }
3916 return false;
3917}
3918
3920{
3921 leftv h=args;
3922 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3923 {
3924 ideal h1= (ideal)h->Data();
3925 res->rtyp =IDEAL_CMD;
3926 res->data =IsSimplex(h1);
3927 }
3928 return false;
3929}
3930
3932{
3933 leftv h=args;
3934 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3935 {
3936 ideal h1= (ideal)h->Data();
3937 h = h->next;
3938 if((h != NULL)&&(h->Typ() == POLY_CMD))
3939 {
3940 poly p= (poly)h->Data();
3941 h = h->next;
3942 if((h != NULL)&&(h->Typ() == INT_CMD))
3943 {
3944 int d= (int)(long)h->Data();
3945 res->rtyp =IDEAL_CMD;
3946 res->data =triangulations1(h1, p, d);
3947 }
3948 }
3949 }
3950 return false;
3951}
3952
3954{
3955 leftv h=args;
3956 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3957 {
3958 ideal h1= (ideal)h->Data();
3959 h = h->next;
3960 if((h != NULL)&&(h->Typ() == POLY_CMD))
3961 {
3962 poly p= (poly)h->Data();
3963 h = h->next;
3964 if((h != NULL)&&(h->Typ() == POLY_CMD))
3965 {
3966 poly q= (poly)h->Data();
3967 h = h->next;
3968 if((h != NULL)&&(h->Typ() == INT_CMD))
3969 {
3970 int d= (int)(long)h->Data();
3971 res->rtyp =IDEAL_CMD;
3972 res->data =triangulations2(h1,p,q,d);
3973 }
3974 }
3975 }
3976 }
3977 return false;
3978}
3979
3981{
3982 leftv h=args;
3983 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
3984 {
3985 ideal h1= (ideal)h->Data();
3986 h = h->next;
3987 if((h != NULL)&&(h->Typ() == POLY_CMD))
3988 {
3989 poly p= (poly)h->Data();
3990 h = h->next;
3991 if((h != NULL)&&(h->Typ() == POLY_CMD))
3992 {
3993 poly q= (poly)h->Data();
3994 h = h->next;
3995 if((h != NULL)&&(h->Typ() == POLY_CMD))
3996 {
3997 poly g= (poly)h->Data();
3998 h = h->next;
3999 if((h != NULL)&&(h->Typ() == INT_CMD))
4000 {
4001 int d= (int)(long)h->Data();
4002 res->rtyp =IDEAL_CMD;
4003 res->data =triangulations3(h1,p,q,g,d);
4004 }
4005 }
4006 }
4007 }
4008 }
4009 return false;
4010}
4011
4013{
4014 leftv h=args;int i;
4015 std::vector<int> bset,bs;
4016 std::vector<std::vector<int> > gset;
4017 if((h != NULL)&&(h->Typ() == INT_CMD))
4018 {
4019 int n= (int)(long)h->Data();
4020 h = h->next;
4021 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4022 {
4023 ideal bi= (ideal)h->Data();
4024 h = h->next;
4025 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4026 {
4027 ideal gi= (ideal)h->Data();
4028 for(i=0;i<IDELEMS(bi);i++)
4029 {
4030 bs=support1(bi->m[i]);
4031 if(bs.size()==1)
4032 bset.push_back(bs[0]);
4033 else if(bs.size()==0)
4034 ;
4035 else
4036 {
4037 WerrorS("Errors in T^1 Equations Solving!");
4038 usleep(1000000);
4039 assert(false);
4040 }
4041
4042 }
4043 gset=supports2(gi);
4044 res->rtyp =INTVEC_CMD;
4045 std::vector<std::vector<int> > vecs=eli2(n,bset,gset);
4046 res->data =Tmat(vecs);
4047 }
4048 }
4049 }
4050 return false;
4051}
4052
4054{
4055 leftv h=args;
4056 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4057 {
4058 ideal h1= (ideal)h->Data();
4059 res->rtyp =IDEAL_CMD;
4060 res->data =trisets(h1);
4061 }
4062 return false;
4063}
4064
4066{
4067 leftv h=args;
4068 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4069 {
4070 ideal h1= (ideal)h->Data();
4071 h = h->next;
4072 if((h != NULL)&&(h->Typ() == POLY_CMD))
4073 {
4074 poly p= (poly)h->Data();
4075 res->rtyp =INT_CMD;
4076 res->data =(void *)(long)valency(h1,p);
4077 }
4078 }
4079 return false;
4080}
4081
4083{
4084 leftv h=args;
4085 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4086 {
4087 ideal h1= (ideal)h->Data();
4088 h = h->next;
4089 if((h != NULL)&&(h->Typ() == POLY_CMD))
4090 {
4091 poly p= (poly)h->Data();
4092 h = h->next;
4093 if((h != NULL)&&(h->Typ() == POLY_CMD))
4094 {
4095 poly q= (poly)h->Data();
4096 res->rtyp =IDEAL_CMD;
4097 std::vector<std::vector<int> > vecs=supports(h1);
4098 std::vector<int> pv=support1(p), qv=support1(q);
4099 res->data =idMaken(Nabv(vecs,pv,qv));
4100 }
4101 }
4102 }
4103 return false;
4104}
4105
4107{
4108 leftv h=args;
4109 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4110 {
4111 ideal h1= (ideal)h->Data();
4112 h = h->next;
4113 if((h != NULL)&&(h->Typ() == POLY_CMD))
4114 {
4115 poly p= (poly)h->Data();
4116 h = h->next;
4117 if((h != NULL)&&(h->Typ() == POLY_CMD))
4118 {
4119 poly q= (poly)h->Data();
4120 res->rtyp =IDEAL_CMD;
4121 std::vector<std::vector<int> > vecs=supports(h1), sbv,tnbr;
4122 std::vector<int> pv=support1(p), qv=support1(q);
4123 std::vector<std::vector<int> > nvs=Nabv(vecs, pv, qv);
4124 ideal sub=psubset(q);
4125 sbv=supports(sub);
4126 std::vector<int> tnv =tnab(vecs,nvs,sbv);
4127 for(unsigned i=0;i<tnv.size();i++)
4128 {
4129 tnbr.push_back(nvs[tnv[i]]);
4130 }
4131 res->data =idMaken(tnbr);
4132 }
4133 }
4134 }
4135 return false;
4136}
4137
4139{
4140 leftv h=args;
4141 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4142 {
4143 ideal h1= (ideal)h->Data();
4144 h = h->next;
4145 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4146 {
4147 ideal h2= (ideal)h->Data();
4148 res->rtyp =INT_CMD;
4149 std::vector<std::vector<int> > vs1=supports(h1), vs2=supports(h2);
4150 res->data =(void *)(long)(vsIntersection(vs1, vs2).size());
4151 }
4152 }
4153 return false;
4154}
4155
4157{
4158 leftv h=args;
4159 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4160 {
4161 ideal h1= (ideal)h->Data();
4162 h = h->next;
4163 if((h != NULL)&&(h->Typ() == POLY_CMD))
4164 {
4165 poly p= (poly)h->Data();
4166 h = h->next;
4167 if((h != NULL)&&(h->Typ() == POLY_CMD))
4168 {
4169 poly q= (poly)h->Data();
4170 res->rtyp =IDEAL_CMD;
4171 res->data =idMaken(Mabv(h1,p,q));
4172 }
4173 }
4174 }
4175 return false;
4176}
4177
4179{
4180 leftv h=args;
4181 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4182 {
4183 ideal h1= (ideal)h->Data();
4184 h = h->next;
4185 if((h != NULL)&&(h->Typ() == POLY_CMD))
4186 {
4187 poly p= (poly)h->Data();
4188 h = h->next;
4189 if((h != NULL)&&(h->Typ() == POLY_CMD))
4190 {
4191 poly q= (poly)h->Data();
4192 std::vector<std::vector<int> > hvs=supports(h1), nv, ntvs;
4193 std::vector<int> av=support1(p), bv=support1(q);
4194 nv=Nabv(hvs,av,bv);
4195 ntvs=nabtv( hvs, nv, av, bv);
4196 std::vector<std::vector<poly> > pvs=idMakei(nv,ntvs);
4197 ideal gens=idInit(1,1);
4198 for(unsigned i=0;i<pvs.size();i++)
4199 {
4200 idInsertPoly(gens,pvs[i][0]);
4201 idInsertPoly(gens,pvs[i][1]);
4202 }
4203 idSkipZeroes(gens);
4204 res->rtyp =IDEAL_CMD;
4205 res->data =gens;
4206 }
4207 }
4208 }
4209 return false;
4210}
4211
4213{
4214 leftv h=args;
4215 if((h != NULL)&&(h->Typ() == POLY_CMD))
4216 {
4217 poly a= (poly)h->Data();
4218 h = h->next;
4219 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4220 {
4221 ideal Xo= (ideal)h->Data();
4222 h = h->next;
4223 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4224 {
4225 ideal Sigma= (ideal)h->Data();
4226 h = h->next;
4227 if((h != NULL)&&(h->Typ() == INT_CMD))
4228 {
4229 int vert= (int)(long)h->Data();
4230 h = h->next;
4231 if((h != NULL)&&(h->Typ() == INT_CMD))
4232 {
4233 int ord= (int)(long)h->Data();
4234 res->rtyp =IDEAL_CMD;
4235 res->data =idMaken(links_new(a, Xo, Sigma, vert, ord));
4236 }
4237 }
4238 }
4239 }
4240 }
4241 return false;
4242}
4243
4245{
4246 leftv h=args;
4247 if((h != NULL)&&(h->Typ() == POLY_CMD))
4248 {
4249 poly p= (poly)h->Data();
4250 h = h->next;
4251 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4252 {
4253 ideal h1= (ideal)h->Data();
4254 res->rtyp =INT_CMD;
4255 res->data =(void *)(long)existIn(p, h1);
4256 }
4257 }
4258 return false;
4259}
4260
4262{
4263 leftv h=args;
4264 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4265 {
4266 ideal h1= (ideal)h->Data();
4267 h = h->next;
4268 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4269 {
4270 ideal h2= (ideal)h->Data();
4271 res->rtyp =IDEAL_CMD;
4272 res->data =idMaken(p_constant(h1,h2));
4273 }
4274 }
4275 return false;
4276}
4277
4279{
4280 leftv h=args;
4281 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4282 {
4283 ideal h1= (ideal)h->Data();
4284 res->rtyp =IDEAL_CMD;
4285 res->data =idMaken(p_change(h1));
4286 }
4287 return false;
4288}
4289
4291{
4292 leftv h=args;
4293 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4294 {
4295 ideal h1= (ideal)h->Data();
4296 h = h->next;
4297 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4298 {
4299 ideal h2= (ideal)h->Data();
4300 res->rtyp =IDEAL_CMD;
4301 res->data =idMaken(p_new(h1,h2));
4302 }
4303 }
4304 return false;
4305}
4306
4308{
4309 leftv h=args;
4310 if((h != NULL)&&(h->Typ() == POLY_CMD))
4311 {
4312 poly p= (poly)h->Data();
4313 res->rtyp =INT_CMD;
4314 res->data =(void *)(long)(support1(p).size());
4315 }
4316 return false;
4317}
4318
4320{
4321 leftv h=args;
4322 if((h != NULL)&&(h->Typ() == POLY_CMD))
4323 {
4324 poly p= (poly)h->Data();
4325 res->rtyp =IDEAL_CMD;
4326 res->data =idMaken(bsubsets_1(p));
4327 }
4328 return false;
4329}
4330
4332{
4333 leftv h=args;
4334 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4335 {
4336 ideal h1= (ideal)h->Data();
4337 h = h->next;
4338 if((h != NULL)&&(h->Typ() == POLY_CMD))
4339 {
4340 poly p= (poly)h->Data();
4341 res->rtyp =IDEAL_CMD;
4342 res->data =idMinusp(h1, p);
4343 }
4344 }
4345 return false;
4346}
4347
4349{
4350 leftv h=args;
4351 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4352 {
4353 ideal h1= (ideal)h->Data();
4354 h = h->next;
4355 if((h != NULL)&&(h->Typ() == POLY_CMD))
4356 {
4357 poly p= (poly)h->Data();
4358 std::vector<std::vector<int> > st=star(p, h1);
4359 std::vector<std::vector<int> > hvs=supports(h1);
4360 std::vector<std::vector<int> > re= vsMinusvs(hvs, st);
4361 res->rtyp =IDEAL_CMD;
4362 res->data =idMaken(re);
4363 }
4364 }
4365 return false;
4366}
4367
4369{
4370 leftv h=args;
4371 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4372 {
4373 ideal h1= (ideal)h->Data();
4374 h = h->next;
4375 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4376 {
4377 ideal h2= (ideal)h->Data();
4378 res->rtyp =IDEAL_CMD;
4379 res->data =c_New(h1, h2);
4380 }
4381 }
4382 return false;
4383}
4384
4386{
4387 leftv h=args;
4388 if((h != NULL)&&(h->Typ() == POLY_CMD))
4389 {
4390 poly p= (poly)h->Data();
4391 h = h->next;
4392 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4393 {
4394 ideal h1= (ideal)h->Data();
4395 res->rtyp =IDEAL_CMD;
4396 res->data =idMaken(star(p, h1));
4397 }
4398 }
4399 return false;
4400}
4401
4403{
4404 leftv h=args;
4405 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4406 {
4407 ideal h2= (ideal)h->Data();
4408 h = h->next;
4409 if((h != NULL)&&(h->Typ() == POLY_CMD))
4410 {
4411 poly p= (poly)h->Data();
4412 res->rtyp =IDEAL_CMD;
4413 res->data =idMaken(stellarsub(p, h2));
4414 }
4415 }
4416 return false;
4417}
4418
4420{
4421 leftv h=args;
4422 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4423 {
4424 ideal h1= (ideal)h->Data();
4425 h = h->next;
4426 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4427 {
4428 ideal h2= (ideal)h->Data();
4429 res->rtyp =IDEAL_CMD;
4430 res->data =idmodulo(h1, h2);
4431 }
4432 }
4433 return false;
4434}
4435
4437{
4438 leftv h=args;
4439 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4440 {
4441 ideal h1= (ideal)h->Data();
4442 h = h->next;
4443 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4444 {
4445 ideal h2= (ideal)h->Data();
4446 res->rtyp =IDEAL_CMD;
4447 res->data =idMinus(h1, h2);
4448 }
4449 }
4450 return false;
4451}
4452
4454{
4455 leftv h=args;
4456 if((h != NULL)&&(h->Typ() == POLY_CMD))
4457 {
4458 poly p= (poly)h->Data();
4459 h = h->next;
4460 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4461 {
4462 ideal h1= (ideal)h->Data();
4463 h = h->next;
4464 if((h != NULL)&&(h->Typ() == POLY_CMD))
4465 {
4466 poly a= (poly)h->Data();
4467 h = h->next;
4468 if((h != NULL)&&(h->Typ() == POLY_CMD))
4469 {
4470 poly b= (poly)h->Data();
4471 res->rtyp =INT_CMD;
4472 res->data =(void *)(long)isoNum(p, h1, a, b);
4473 }
4474 }
4475 }
4476 }
4477 return false;
4478}
4479
4481{
4482 leftv h=args;
4483 if((h != NULL)&&(h->Typ() == POLY_CMD))
4484 {
4485 poly p= (poly)h->Data();
4486 h = h->next;
4487 if((h != NULL)&&(h->Typ() == POLY_CMD))
4488 {
4489 poly q= (poly)h->Data();
4490 h = h->next;
4491 if((h != NULL)&&(h->Typ() == POLY_CMD))
4492 {
4493 poly f= (poly)h->Data();
4494 h = h->next;
4495 if((h != NULL)&&(h->Typ() == POLY_CMD))
4496 {
4497 poly g= (poly)h->Data();
4498 h = h->next;
4499 if((h != NULL)&&(h->Typ() == POLY_CMD))
4500 {
4501 poly a= (poly)h->Data();
4502 h = h->next;
4503 if((h != NULL)&&(h->Typ() == POLY_CMD))
4504 {
4505 poly b= (poly)h->Data();
4506 res->rtyp =INT_CMD;
4507 res->data =(void *)(long)ifIso(p,q,f,g, a, b);
4508 }
4509 }
4510 }
4511 }
4512 }
4513 }
4514 return false;
4515}
4516
4518{
4519 leftv h=args;
4520 if((h != NULL)&&(h->Typ() == POLY_CMD))
4521 {
4522 poly p= (poly)h->Data();
4523 h = h->next;
4524 if((h != NULL)&&(h->Typ() == INT_CMD))
4525 {
4526 int num= (int)(long)h->Data();
4527 res->rtyp =INT_CMD;
4528 res->data =(void *)(long)redefinedeg( p, num);
4529 }
4530 }
4531 return false;
4532}
4533
4535{
4536 leftv h=args;
4537 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4538 {
4539 ideal h1= (ideal)h->Data();
4540 res->rtyp =IDEAL_CMD;
4541 res->data =complementsimplex(h1);
4542 }
4543 return false;
4544}
4545
4547{
4548 leftv h=args;
4549 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4550 {
4551 ideal h1= (ideal)h->Data();
4552 res->rtyp =INT_CMD;
4553 res->data =(void *)(long)dim_sim(h1);
4554 }
4555 return false;
4556}
4557
4559{
4560 leftv h=args;
4561 if((h != NULL)&&(h->Typ() == IDEAL_CMD))
4562 {
4563 ideal h1= (ideal)h->Data();
4564 h = h->next;
4565 if((h != NULL)&&(h->Typ() == INT_CMD))
4566 {
4567 int num= (int)(long)h->Data();
4568 res->rtyp =INT_CMD;
4569 res->data =(void *)(long)num4dim( h1, num);
4570 }
4571 }
4572 return false;
4573}
4574
4575/**************************************interface T2****************************************/
4576
4578{
4579 p->iiAddCproc("","mg",FALSE,idsr);
4580 p->iiAddCproc("","gd",FALSE,gd);
4581 p->iiAddCproc("","findbset",FALSE,fb);
4582 p->iiAddCproc("","findaset",FALSE,fa);
4583 p->iiAddCproc("","fgp",FALSE,fgp);
4584 p->iiAddCproc("","fgpl",FALSE,fgpl);
4585 p->iiAddCproc("","idcomplement",FALSE,idcomplement);
4586 p->iiAddCproc("","genst",FALSE,genstt);
4587 p->iiAddCproc("","sgp",FALSE,sgp);
4588 p->iiAddCproc("","sgpl",FALSE,sgpl);
4589 p->iiAddCproc("","Links",FALSE,Links);
4590 p->iiAddCproc("","eqsolve1",FALSE,eqsolve1);
4591 p->iiAddCproc("","pb",FALSE,pb);
4592 p->iiAddCproc("","pa",FALSE,pa);
4593 p->iiAddCproc("","makeSimplex",FALSE,makeSimplex);
4594 p->iiAddCproc("","isSim",FALSE,isSim);
4595 p->iiAddCproc("","nfaces1",FALSE,nfaces1);
4596 p->iiAddCproc("","nfaces2",FALSE,nfaces2);
4597 p->iiAddCproc("","nfaces3",FALSE,nfaces3);
4598 p->iiAddCproc("","comedg",FALSE,comedg);
4599 p->iiAddCproc("","tsets",FALSE,tsets);
4600 p->iiAddCproc("","valency",FALSE,Valency);
4601 p->iiAddCproc("","nab",FALSE,nabvl);
4602 p->iiAddCproc("","tnab",FALSE,tnabvl);
4603 p->iiAddCproc("","mab",FALSE,mabvl);
4604 p->iiAddCproc("","SRideal",FALSE,SRideal);
4605 p->iiAddCproc("","Linkn",FALSE,linkn);
4606 p->iiAddCproc("","Existb",FALSE,existsub);
4607 p->iiAddCproc("","pConstant",FALSE,pConstant);
4608 p->iiAddCproc("","pChange",FALSE,pChange);
4609 p->iiAddCproc("","pNew",FALSE,p_New);
4610 p->iiAddCproc("","pSupport",FALSE,support);
4611 p->iiAddCproc("","psMinusp",FALSE,psMinusp);
4612 p->iiAddCproc("","cNew",FALSE,cNew);
4613 p->iiAddCproc("","isoNumber",FALSE,isoNumber);
4614 p->iiAddCproc("","vsInsec",FALSE,vsIntersec);
4615 p->iiAddCproc("","getnabt",FALSE,nabtvl);
4616 p->iiAddCproc("","idmodulo",FALSE,idModulo);
4617 p->iiAddCproc("","ndegree",FALSE,newDegree);
4618 p->iiAddCproc("","nonf2f",FALSE,nonf2f);
4619 p->iiAddCproc("","ifIsom",FALSE,ifIsomorphism);
4620 p->iiAddCproc("","stellarsubdivision",FALSE,stellarsubdivision);
4621 p->iiAddCproc("","star",FALSE,stars);
4622 p->iiAddCproc("","numdim",FALSE,numdim);
4623 p->iiAddCproc("","dimsim",FALSE,dimsim);
4624 p->iiAddCproc("","bprime",FALSE,bprime);
4625 p->iiAddCproc("","remainpart",FALSE,stellarremain);
4626 p->iiAddCproc("","idminus",FALSE,idminus);
4627 p->iiAddCproc("","time1",FALSE,t1h);
4628}
4629
4630extern "C" int SI_MOD_INIT(cohomo)(SModulFunctions* p)
4631{
4633 return MAX_TOK;
4634}
4635#endif
4636#endif
int BOOLEAN
Definition: auxiliary.h:87
#define FALSE
Definition: auxiliary.h:96
int size(const CanonicalForm &f, const Variable &v)
int size ( const CanonicalForm & f, const Variable & v )
Definition: cf_ops.cc:600
CanonicalForm FACTORY_PUBLIC pp(const CanonicalForm &)
CanonicalForm pp ( const CanonicalForm & f )
Definition: cf_gcd.cc:676
CanonicalForm num(const CanonicalForm &f)
Variable mvar(const CanonicalForm &f)
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:56
int l
Definition: cfEzgcd.cc:100
int m
Definition: cfEzgcd.cc:128
int i
Definition: cfEzgcd.cc:132
int p
Definition: cfModGcd.cc:4078
g
Definition: cfModGcd.cc:4090
CanonicalForm cf
Definition: cfModGcd.cc:4083
CanonicalForm gp
Definition: cfModGcd.cc:4102
CanonicalForm b
Definition: cfModGcd.cc:4103
bool solve(int **extmat, int nrows, int ncols)
Definition: cf_linsys.cc:504
FILE * f
Definition: checklibs.c:9
Definition: idrec.h:35
Definition: intvec.h:23
Class used for (list of) interpreter objects.
Definition: subexpr.h:83
void * CopyD(int t)
Definition: subexpr.cc:710
int rtyp
Definition: subexpr.h:91
void * Data()
Definition: subexpr.cc:1162
void Init()
Definition: subexpr.h:107
void * data
Definition: subexpr.h:88
Definition: lists.h:24
sleftv * m
Definition: lists.h:46
Coefficient rings, fields and other domains suitable for Singular polynomials.
static FORCE_INLINE long n_Int(number &n, const coeffs r)
conversion of n to an int; 0 if not possible in Z/pZ: the representing int lying in (-p/2 ....
Definition: coeffs.h:544
static FORCE_INLINE coeffs nCopyCoeff(const coeffs r)
"copy" coeffs, i.e. increment ref
Definition: coeffs.h:430
static BOOLEAN fa(leftv res, leftv args)
Definition: cohomo.cc:3764
static std::vector< int > vecbase1(int num, std::vector< int > oset)
Definition: cohomo.cc:1169
static BOOLEAN pa(leftv res, leftv args)
Definition: cohomo.cc:3723
static BOOLEAN tsets(leftv res, leftv args)
Definition: cohomo.cc:4053
static ideal T_1h(ideal h)
Definition: cohomo.cc:3549
static std::vector< int > fvarsvalue(int vnum, std::vector< int > fvars)
Definition: cohomo.cc:1108
static std::vector< int > commonedge(poly p, poly q)
Definition: cohomo.cc:3042
static std::vector< int > v_minus(std::vector< int > v1, std::vector< int > v2)
Definition: cohomo.cc:3446
static bool IsinL(int a, std::vector< int > vec)
Definition: cohomo.cc:132
static BOOLEAN tnabvl(leftv res, leftv args)
Definition: cohomo.cc:4106
BOOLEAN nfaces1(leftv res, leftv args)
Definition: cohomo.cc:3931
static std::vector< int > subspace1(std::vector< std::vector< int > > mv, std::vector< int > bv)
Definition: cohomo.cc:1640
VAR clock_t t_construct
Definition: cohomo.cc:2692
static std::vector< int > tnab(std::vector< std::vector< int > > hvs, std::vector< std::vector< int > > nvs, std::vector< std::vector< int > > bvs)
Definition: cohomo.cc:2172
static BOOLEAN cNew(leftv res, leftv args)
Definition: cohomo.cc:4368
static std::vector< std::vector< int > > listsinsertlist(std::vector< std::vector< int > > gset, int a, int b)
Definition: cohomo.cc:1561
static std::vector< std::vector< int > > vsMinusvs(std::vector< std::vector< int > > vs1, std::vector< std::vector< int > > vs2)
Definition: cohomo.cc:3226
static BOOLEAN comedg(leftv res, leftv args)
Definition: cohomo.cc:3694
static poly pMaken(std::vector< int > vbase)
Definition: cohomo.cc:465
static std::vector< std::vector< int > > value1l(std::vector< std::vector< int > > mvs, std::vector< std::vector< int > > lks, std::vector< std::vector< int > > vecs, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2659
VAR clock_t t_start
Definition: cohomo.cc:2692
static BOOLEAN idModulo(leftv res, leftv args)
Definition: cohomo.cc:4419
static ideal idMaken(std::vector< std::vector< int > > vecs)
Definition: cohomo.cc:479
static std::vector< int > vecIntersection(std::vector< int > p, std::vector< int > q)
Definition: cohomo.cc:147
static std::vector< std::vector< int > > penface(poly p, poly q, poly g, int vert)
Definition: cohomo.cc:3109
static ideal idMake(std::vector< std::vector< int > > vecs)
Definition: cohomo.cc:358
static std::vector< std::vector< int > > vsMinusv(std::vector< std::vector< int > > vecs, std::vector< int > vec)
Definition: cohomo.cc:225
static intvec * gradedpiece1n(ideal h, poly a, poly b)
Definition: cohomo.cc:2313
static std::vector< std::vector< int > > value2l(std::vector< std::vector< int > > mvs, std::vector< std::vector< int > > lks, std::vector< std::vector< int > > mts, std::vector< std::vector< int > > lkts, std::vector< std::vector< int > > vecs, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2792
static BOOLEAN isoNumber(leftv res, leftv args)
Definition: cohomo.cc:4453
static BOOLEAN makeSimplex(leftv res, leftv args)
Definition: cohomo.cc:3735
static bool vInvsl(std::vector< int > vec, std::vector< std::vector< int > > vecs)
Definition: cohomo.cc:186
static std::vector< std::vector< int > > mabtv(std::vector< std::vector< int > > hvs, std::vector< std::vector< int > > Mv, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:1983
static BOOLEAN fgp(leftv res, leftv args)
Definition: cohomo.cc:3786
static bool tNab(std::vector< std::vector< int > > hvs, std::vector< int > pv, std::vector< std::vector< int > > bvs)
Definition: cohomo.cc:2157
static int valency(ideal h, poly p)
Definition: cohomo.cc:3178
static bool condition3for2(std::vector< std::vector< int > > hvs, std::vector< int > pv, std::vector< int > qv, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:1777
static void lpsprint(std::vector< std::vector< poly > > pvs)
Definition: cohomo.cc:112
static std::vector< std::vector< int > > value2(std::vector< std::vector< int > > mvs, std::vector< std::vector< int > > nvs, std::vector< std::vector< int > > mts, std::vector< std::vector< int > > nts, std::vector< std::vector< int > > vecs, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2459
static std::vector< int > phimagel(std::vector< int > fv, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2649
static BOOLEAN stars(leftv res, leftv args)
Definition: cohomo.cc:4385
static std::vector< int > keeporder(std::vector< int > vec)
Definition: cohomo.cc:1039
static std::vector< std::vector< int > > p_new(ideal Xo, ideal Sigma)
Definition: cohomo.cc:3262
static BOOLEAN newDegree(leftv res, leftv args)
Definition: cohomo.cc:4517
static void T1(ideal h)
Definition: cohomo.cc:2380
static bool condition2for2(std::vector< std::vector< int > > hvs, std::vector< int > pv, std::vector< int > qv, std::vector< int > sv, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:1764
static std::vector< std::vector< int > > Nabv(std::vector< std::vector< int > > hvs, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2106
static ideal trisets(ideal h)
Definition: cohomo.cc:2969
static int pvert(poly p)
Definition: cohomo.cc:541
static int idvert(ideal h)
Definition: cohomo.cc:522
static void firstorderdef_setup(SModulFunctions *p)
Definition: cohomo.cc:4577
static ideal SimFacset(poly p)
Definition: cohomo.cc:787
static std::vector< std::vector< int > > links_new(poly a, ideal Xo, ideal Sigma, int vert, int ord)
Definition: cohomo.cc:3361
static std::vector< std::vector< poly > > idMakei(std::vector< std::vector< int > > mv, std::vector< std::vector< int > > vecs)
Definition: cohomo.cc:1671
static ideal psubset(poly p)
Definition: cohomo.cc:1539
static BOOLEAN bprime(leftv res, leftv args)
Definition: cohomo.cc:4319
static void listprint(std::vector< int > vec)
Definition: cohomo.cc:49
static BOOLEAN pConstant(leftv res, leftv args)
Definition: cohomo.cc:4261
static bool vInp(int m, poly p)
Definition: cohomo.cc:405
static BOOLEAN dimsim(leftv res, leftv args)
Definition: cohomo.cc:4546
static std::vector< std::vector< int > > ofindbases(int num, std::vector< int > bset, std::vector< std::vector< int > > gset)
Definition: cohomo.cc:1227
static BOOLEAN fgpl(leftv res, leftv args)
Definition: cohomo.cc:3808
static bool p_Ifsfree(poly P)
Definition: cohomo.cc:643
static int redefinedeg(poly p, int num)
Definition: cohomo.cc:1315
static intvec * edgemat(poly p, poly q)
Definition: cohomo.cc:3053
static BOOLEAN psMinusp(leftv res, leftv args)
Definition: cohomo.cc:4331
static bool nabtconditionv(std::vector< std::vector< int > > hvs, std::vector< int > pv, std::vector< int > qv)
Definition: cohomo.cc:2123
static BOOLEAN SRideal(leftv res, leftv args)
Definition: cohomo.cc:3600
static std::vector< std::vector< int > > phi1(poly a, ideal Sigma)
Definition: cohomo.cc:3328
static bool condition2for2nv(std::vector< std::vector< int > > hvs, std::vector< int > pv, std::vector< int > qv, std::vector< int > fv)
Definition: cohomo.cc:2410
static std::vector< int > freevars(int n, std::vector< int > bset, std::vector< std::vector< int > > gset)
Definition: cohomo.cc:1085
static BOOLEAN stellarsubdivision(leftv res, leftv args)
Definition: cohomo.cc:4402
static BOOLEAN idcomplement(leftv res, leftv args)
Definition: cohomo.cc:3612
static std::vector< std::vector< int > > phi2(poly a, ideal Xo, ideal Sigma)
Definition: cohomo.cc:3344
VAR clock_t t_begin
Definition: cohomo.cc:2692
static int isoNum(poly p, ideal I, poly a, poly b)
Definition: cohomo.cc:3401
static bool condition1for2(std::vector< int > pv, std::vector< int > qv, std::vector< int > bv)
Definition: cohomo.cc:1751
static BOOLEAN isSim(leftv res, leftv args)
Definition: cohomo.cc:3919
static ideal IsSimplex(ideal h)
Definition: cohomo.cc:832
static BOOLEAN gd(leftv res, leftv args)
Definition: cohomo.cc:3677
static std::vector< std::vector< int > > vsUnion(std::vector< std::vector< int > > vs1, std::vector< std::vector< int > > vs2)
Definition: cohomo.cc:239
static BOOLEAN nabvl(leftv res, leftv args)
Definition: cohomo.cc:4082
static BOOLEAN fb(leftv res, leftv args)
Definition: cohomo.cc:3711
static ideal idmodulo(ideal h1, ideal h2)
Definition: cohomo.cc:375
VAR clock_t t_mark
Definition: cohomo.cc:2692
static ideal idadda(ideal h1, ideal h2)
Definition: cohomo.cc:808
static intvec * gradedpiece2n(ideal h, poly a, poly b)
Definition: cohomo.cc:2528
static ideal idMake3(std::vector< std::vector< int > > vecs)
Definition: cohomo.cc:1601
static ideal complementsimplex(ideal h)
Definition: cohomo.cc:853
static ideal c_New(ideal Io, ideal sig)
Definition: cohomo.cc:3293
static void T2(ideal h)
Definition: cohomo.cc:2607
static std::vector< int > make1(int n)
Definition: cohomo.cc:1198
static ideal qringadd(ideal h1, ideal h2, int deg)
Definition: cohomo.cc:735
static BOOLEAN eqsolve1(leftv res, leftv args)
Definition: cohomo.cc:4012
static std::vector< std::vector< int > > vAbsorb(std::vector< int > bset, std::vector< std::vector< int > > gset)
Definition: cohomo.cc:1125
static ideal finda(ideal h, poly S, int ddeg)
Definition: cohomo.cc:934
static std::vector< std::vector< int > > soleli1(std::vector< std::vector< int > > eqs)
Definition: cohomo.cc:1052
static std::vector< int > vMake(poly p)
Definition: cohomo.cc:420
static std::vector< int > subspacet1(int num, std::vector< std::vector< int > > ntvs)
Definition: cohomo.cc:1942
static std::vector< std::vector< int > > vsMake(ideal h)
Definition: cohomo.cc:439
static BOOLEAN numdim(leftv res, leftv args)
Definition: cohomo.cc:4558
static bool vEvl(std::vector< int > vec1, std::vector< int > vec2)
Definition: cohomo.cc:174
static BOOLEAN vsIntersec(leftv res, leftv args)
Definition: cohomo.cc:4138
static BOOLEAN support(leftv res, leftv args)
Definition: cohomo.cc:4307
static void id_print(ideal h)
Definition: cohomo.cc:84
static std::vector< std::vector< int > > bsubsets_1(poly b)
Definition: cohomo.cc:3532
static BOOLEAN genstt(leftv res, leftv args)
Definition: cohomo.cc:3835
static void TimeShow(clock_t t_construct, clock_t t_solve, clock_t t_value, clock_t t_total)
Definition: cohomo.cc:2695
static intvec * dmat(poly a, poly b)
Definition: cohomo.cc:3659
static BOOLEAN nonf2f(leftv res, leftv args)
Definition: cohomo.cc:4534
static std::vector< std::vector< int > > stellarsub(poly a, ideal h)
Definition: cohomo.cc:3499
static void equmab(int num)
Definition: cohomo.cc:1618
static BOOLEAN sgpl(leftv res, leftv args)
Definition: cohomo.cc:3879
static std::vector< int > support1(poly p)
Definition: cohomo.cc:271
static bool nabconditionv(std::vector< std::vector< int > > hvs, std::vector< int > pv, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2092
static std::vector< int > support2(poly p)
Definition: cohomo.cc:304
static std::vector< std::vector< int > > p_change(ideal Sigma)
Definition: cohomo.cc:3255
static std::vector< std::vector< int > > vecqring(std::vector< std::vector< int > > vec1, std::vector< std::vector< int > > vec2)
Definition: cohomo.cc:454
static std::vector< int > gdegree(poly a, poly b)
Definition: cohomo.cc:3456
static ideal id_sfmon(ideal h)
Definition: cohomo.cc:679
static ideal triangulations1(ideal h, poly p, int vert)
Definition: cohomo.cc:3004
VAR clock_t t_value
Definition: cohomo.cc:2692
static int id_maxdeg(ideal h)
Definition: cohomo.cc:744
static std::vector< int > vertset(std::vector< std::vector< int > > vecs)
Definition: cohomo.cc:1404
static ideal id_complement(ideal h)
Definition: cohomo.cc:698
static ideal p_a(ideal h)
Definition: cohomo.cc:1335
static std::vector< std::vector< int > > tetraface(poly p, poly q, int vert)
Definition: cohomo.cc:3072
static std::vector< std::vector< int > > eli2(int num, std::vector< int > bset, std::vector< std::vector< int > > gset)
Definition: cohomo.cc:1254
static ideal triangulations3(ideal h, poly p, poly q, poly g, int vert)
Definition: cohomo.cc:3153
static std::vector< std::vector< int > > vs_subsets(std::vector< std::vector< int > > vs)
Definition: cohomo.cc:3236
static std::vector< int > vecUnion(std::vector< int > vec1, std::vector< int > vec2)
Definition: cohomo.cc:200
static BOOLEAN sgp(leftv res, leftv args)
Definition: cohomo.cc:3857
static intvec * Tmat(std::vector< std::vector< int > > vecs)
Definition: cohomo.cc:2233
static std::vector< std::vector< int > > canonicalbase(int n)
Definition: cohomo.cc:1843
static BOOLEAN idsr(leftv res, leftv args)
Definition: cohomo.cc:3637
static std::vector< std::vector< int > > minisolve(std::vector< std::vector< int > > solve, std::vector< int > index)
Definition: cohomo.cc:2291
static ideal genst(ideal h, poly a, poly b)
Definition: cohomo.cc:2511
static BOOLEAN linkn(leftv res, leftv args)
Definition: cohomo.cc:4212
static std::vector< int > vecMinus(std::vector< int > vec1, std::vector< int > vec2)
Definition: cohomo.cc:212
static std::vector< std::vector< int > > Mabv(ideal h, poly a, poly b)
Definition: cohomo.cc:975
static std::vector< int > gensindex(ideal M, ideal ids)
Definition: cohomo.cc:2260
static BOOLEAN Valency(leftv res, leftv args)
Definition: cohomo.cc:4065
static std::vector< int > makeequation(int i, int j, int t)
Definition: cohomo.cc:1570
static bool mabconditionv(std::vector< std::vector< int > > hvs, std::vector< int > pv, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:963
static std::vector< std::vector< int > > boundary(poly a)
Definition: cohomo.cc:3489
static ideal idMinus(ideal h1, ideal h2)
Definition: cohomo.cc:617
int SI_MOD_INIT() cohomo(SModulFunctions *p)
Definition: cohomo.cc:4630
static int dim_sim(ideal h)
Definition: cohomo.cc:876
static std::vector< std::vector< int > > star(poly a, ideal h)
Definition: cohomo.cc:3473
static std::vector< int > findalpha(std::vector< std::vector< int > > mv, std::vector< int > bv)
Definition: cohomo.cc:1927
static BOOLEAN stellarremain(leftv res, leftv args)
Definition: cohomo.cc:4348
static void lpprint(std::vector< poly > pv)
Definition: cohomo.cc:97
static ideal getpresolve(ideal h)
Definition: cohomo.cc:1796
static BOOLEAN nabtvl(leftv res, leftv args)
Definition: cohomo.cc:4178
static std::vector< int > make0(int n)
Definition: cohomo.cc:1185
static std::vector< std::vector< int > > subspacetn(std::vector< std::vector< int > > N, std::vector< int > tN, std::vector< std::vector< int > > ntvs)
Definition: cohomo.cc:2440
static int pcoef(poly p, int m)
Definition: cohomo.cc:386
static ideal triangulations2(ideal h, poly p, poly q, int vert)
Definition: cohomo.cc:3094
static std::vector< std::vector< int > > getvector(ideal h, int n)
Definition: cohomo.cc:1863
VAR clock_t t_total
Definition: cohomo.cc:2692
static std::vector< std::vector< int > > gpl(ideal h, poly a, poly b)
Definition: cohomo.cc:2704
static BOOLEAN nfaces3(leftv res, leftv args)
Definition: cohomo.cc:3980
static std::vector< std::vector< int > > nabtv(std::vector< std::vector< int > > hvs, std::vector< std::vector< int > > Nv, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2135
static void gradedpiece1(ideal h, poly a, poly b)
Definition: cohomo.cc:1691
static BOOLEAN p_New(leftv res, leftv args)
Definition: cohomo.cc:4290
static BOOLEAN Links(leftv res, leftv args)
Definition: cohomo.cc:3901
static poly pMake(std::vector< int > vbase)
Definition: cohomo.cc:343
static BOOLEAN nfaces2(leftv res, leftv args)
Definition: cohomo.cc:3953
static std::vector< std::vector< int > > value1(std::vector< std::vector< int > > mvs, std::vector< std::vector< int > > nvs, std::vector< std::vector< int > > vecs, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2196
static bool IsInX(poly p, ideal X)
Definition: cohomo.cc:719
static intvec * gradedpiece1nl(ideal h, poly a, poly b, int set)
Definition: cohomo.cc:2768
static std::vector< std::vector< int > > links(poly a, ideal h)
Definition: cohomo.cc:1293
static std::vector< std::vector< int > > supports(ideal h)
Definition: cohomo.cc:287
static BOOLEAN pb(leftv res, leftv args)
Definition: cohomo.cc:3747
VAR clock_t t_solve
Definition: cohomo.cc:2692
static std::vector< int > ofindbases1(int num, int vnum, std::vector< int > bset, std::vector< std::vector< int > > gset)
Definition: cohomo.cc:1211
static std::vector< std::vector< int > > vsIntersection(std::vector< std::vector< int > > vs1, std::vector< std::vector< int > > vs2)
Definition: cohomo.cc:253
static ideal sfreemon(ideal h, int deg)
Definition: cohomo.cc:658
static int existIn(poly b, ideal Xs)
Definition: cohomo.cc:3386
static std::vector< std::vector< int > > supports2(ideal h)
Definition: cohomo.cc:326
static BOOLEAN pChange(leftv res, leftv args)
Definition: cohomo.cc:4278
static std::vector< std::vector< int > > gpl2(ideal h, poly a, poly b)
Definition: cohomo.cc:2845
static std::vector< int > numfree(ideal h)
Definition: cohomo.cc:1824
static std::vector< int > eli1(std::vector< int > eq1, std::vector< int > eq2)
Definition: cohomo.cc:996
static ideal idsrRing(ideal h)
Definition: cohomo.cc:758
static ideal idMinusp(ideal I, poly p)
Definition: cohomo.cc:3429
static std::vector< std::vector< int > > triface(poly p, int vert)
Definition: cohomo.cc:2985
static BOOLEAN idminus(leftv res, leftv args)
Definition: cohomo.cc:4436
static ideal p_b(ideal h, poly a)
Definition: cohomo.cc:1427
static void listsprint(std::vector< std::vector< int > > posMat)
Definition: cohomo.cc:65
static BOOLEAN mabvl(leftv res, leftv args)
Definition: cohomo.cc:4156
static int num4dim(ideal h, int n)
Definition: cohomo.cc:889
static std::vector< int > findalphan(std::vector< std::vector< int > > N, std::vector< int > tN)
Definition: cohomo.cc:2424
static BOOLEAN ifIsomorphism(leftv res, leftv args)
Definition: cohomo.cc:4480
static int ifIso(poly p, poly q, poly f, poly g, poly a, poly b)
Definition: cohomo.cc:3418
static void gradedpiece2(ideal h, poly a, poly b)
Definition: cohomo.cc:2006
static intvec * gradedpiece2nl(ideal h, poly a, poly b)
Definition: cohomo.cc:2913
static std::vector< std::vector< int > > b_subsets(std::vector< int > vec)
Definition: cohomo.cc:496
static std::vector< std::vector< int > > p_constant(ideal Xo, ideal Sigma)
Definition: cohomo.cc:3247
static BOOLEAN existsub(leftv res, leftv args)
Definition: cohomo.cc:4244
static ideal mingens(ideal h, poly a, poly b)
Definition: cohomo.cc:2277
static std::vector< int > phimage(std::vector< int > pv, std::vector< int > av, std::vector< int > bv)
Definition: cohomo.cc:2187
static bool vsubset(std::vector< int > vec1, std::vector< int > vec2)
Definition: cohomo.cc:160
static ideal findb(ideal h)
Definition: cohomo.cc:908
static BOOLEAN t1h(leftv res, leftv args)
Definition: cohomo.cc:3625
static std::vector< poly > pMakei(std::vector< std::vector< int > > mv, std::vector< int > vbase)
Definition: cohomo.cc:1657
static std::vector< std::vector< int > > subspacet(std::vector< std::vector< int > > mv, std::vector< int > bv, std::vector< std::vector< int > > ntvs)
Definition: cohomo.cc:1969
static poly pMake3(std::vector< int > vbase)
Definition: cohomo.cc:1583
#define Print
Definition: emacs.cc:80
const CanonicalForm int s
Definition: facAbsFact.cc:51
CanonicalForm res
Definition: facAbsFact.cc:60
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:39
bool bad
Definition: facFactorize.cc:64
fq_nmod_poly_t * vec
Definition: facHensel.cc:108
int j
Definition: facHensel.cc:110
static int max(int a, int b)
Definition: fast_mult.cc:264
void WerrorS(const char *s)
Definition: feFopen.cc:24
#define VAR
Definition: globaldefs.h:5
@ IDEAL_CMD
Definition: grammar.cc:284
@ POLY_CMD
Definition: grammar.cc:289
@ RING_CMD
Definition: grammar.cc:281
ideal scKBase(int deg, ideal s, ideal Q, intvec *mv)
Definition: hdegree.cc:1449
BOOLEAN idInsertPoly(ideal h1, poly h2)
insert h2 into h1 (if h2 is not the zero polynomial) return TRUE iff h2 was indeed inserted
BOOLEAN idIs0(ideal h)
returns true if h is the zero ideal
ideal idCopy(ideal A)
Definition: ideals.h:60
ideal idAdd(ideal h1, ideal h2)
h1 + h2
Definition: ideals.h:68
#define IMATELEM(M, I, J)
Definition: intvec.h:85
idhdl ggetid(const char *n)
Definition: ipid.cc:581
idhdl enterid(const char *s, int lev, int t, idhdl *root, BOOLEAN init, BOOLEAN search)
Definition: ipid.cc:279
#define IDROOT
Definition: ipid.h:19
#define IDRING(a)
Definition: ipid.h:127
BOOLEAN iiMake_proc(idhdl pn, package pack, leftv args)
Definition: iplib.cc:504
INST_VAR sleftv iiRETURNEXPR
Definition: iplib.cc:474
void rSetHdl(idhdl h)
Definition: ipshell.cc:5126
STATIC_VAR Poly * h
Definition: janet.cc:971
poly kNF(ideal F, ideal Q, poly p, int syzComp, int lazyReduce)
Definition: kstd1.cc:3182
ideal kStd(ideal F, ideal Q, tHomog h, intvec **w, intvec *hilb, int syzComp, int newIdeal, intvec *vw, s_poly_proc_t sp)
Definition: kstd1.cc:2447
void rem(unsigned long *a, unsigned long *q, unsigned long p, int &dega, int degq)
Definition: minpoly.cc:572
#define pNext(p)
Definition: monomials.h:36
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
slists * lists
Definition: mpr_numeric.h:146
The main handler for Singular numbers which are suitable for Singular polynomials.
#define nInit(i)
Definition: numbers.h:24
#define omStrDup(s)
Definition: omAllocDecl.h:263
#define omAlloc(size)
Definition: omAllocDecl.h:210
#define omalloc(size)
Definition: omAllocDecl.h:228
#define NULL
Definition: omList.c:12
static int index(p_Length length, p_Ord ord)
Definition: p_Procs_Impl.h:592
poly p_Subst(poly p, int n, poly e, const ring r)
Definition: p_polys.cc:3958
BOOLEAN p_EqualPolys(poly p1, poly p2, const ring r)
Definition: p_polys.cc:4512
static poly pp_Mult_mm(poly p, poly m, const ring r)
Definition: p_polys.h:1029
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
static BOOLEAN p_DivisibleBy(poly a, poly b, const ring r)
Definition: p_polys.h:1898
void rChangeCurrRing(ring r)
Definition: polys.cc:15
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:13
#define pAdd(p, q)
Definition: polys.h:203
static long pTotaldegree(poly p)
Definition: polys.h:282
#define pSetm(p)
Definition: polys.h:271
#define pSetCoeff(p, n)
deletes old coeff before setting the new one
Definition: polys.h:31
void pWrite(poly p)
Definition: polys.h:308
#define pGetExp(p, i)
Exponent.
Definition: polys.h:41
#define pEqualPolys(p1, p2)
Definition: polys.h:399
#define pSetExp(p, i, v)
Definition: polys.h:42
#define pCopy(p)
return a copy of the poly
Definition: polys.h:185
#define pOne()
Definition: polys.h:315
void PrintS(const char *s)
Definition: reporter.cc:284
void PrintLn()
Definition: reporter.cc:310
ring rDefault(const coeffs cf, int N, char **n, int ord_size, rRingOrder_t *ord, int *block0, int *block1, int **wvhdl, unsigned long bitmask)
Definition: ring.cc:102
ring rCopy(ring r)
Definition: ring.cc:1731
@ ringorder_lp
Definition: ring.h:77
static short rVar(const ring r)
#define rVar(r) (r->N)
Definition: ring.h:592
ideal id_Add(ideal h1, ideal h2, const ring r)
h1 + h2
ideal idInit(int idsize, int rank)
initialise an ideal / module
Definition: simpleideals.cc:35
ideal id_MaxIdeal(const ring r)
initialise the maximal ideal (at 0)
Definition: simpleideals.cc:98
void idSkipZeroes(ideal ide)
gives an ideal/module the minimal possible size
#define IDELEMS(i)
Definition: simpleideals.h:23
#define R
Definition: sirandom.c:27
#define M
Definition: sirandom.c:25
@ testHomog
Definition: structs.h:38
#define assert(A)
Definition: svd_si.h:3
@ INTVEC_CMD
Definition: tok.h:101
@ INT_CMD
Definition: tok.h:96
@ MAX_TOK
Definition: tok.h:218
int dim(ideal I, ring r)
int tdeg(poly p)
Definition: walkSupport.cc:35