source: git/kernel/test.cc @ f323dd1

spielwiese
Last change on this file since f323dd1 was f323dd1, checked in by Hans Schoenemann <hannes@…>, 12 years ago
chg: move feResource to findexec lib
  • Property mode set to 100644
File size: 12.6 KB
Line 
1#include "config.h"
2#include "mod2.h"
3
4#include <omalloc/omalloc.h>
5#include <misc/auxiliary.h>
6#ifdef HAVE_FACTORY
7#include <factory/factory.h> // :(
8#endif
9
10#include <misc/intvec.h>
11#include <misc/int64vec.h>
12#include <misc/mylimits.h>
13#include <misc/options.h>
14
15#include <reporter/reporter.h>
16
17#include <resources/feFopen.h>
18#include <findexec/feResource.h>
19
20#include <coeffs/coeffs.h>
21
22#include <coeffs/si_gmp.h>
23
24#include <polys/kbuckets.h>
25#include <polys/matpol.h>
26#include <polys/mod_raw.h>
27#include <polys/prCopy.h>
28#include <polys/sbuckets.h>
29#include <polys/simpleideals.h>
30#include <polys/weight.h>
31
32#include <polys/monomials/maps.h>
33#include <polys/monomials/monomials.h>
34#include <polys/monomials/p_polys.h>
35#include <polys/monomials/ring.h>
36
37#include <polys/nc/nc.h>
38#include <polys/nc/ncSACache.h>
39#include <polys/nc/ncSAFormula.h>
40#include <polys/nc/ncSAMult.h>
41#include <polys/nc/sca.h>
42#include <polys/nc/summator.h>
43
44
45#include <polys/templates/p_MemAdd.h>
46#include <polys/templates/p_Procs.h>
47
48#include <polys/operations/pShallowCopyDelete.h>
49
50#include <polys/clapsing.h>
51
52
53#ifdef HAVE_FACTORY
54// The following are needed due to FACTORY (e.g. initCanonicalForm)
55int initializeGMP(){ return 1; }
56int mmInit(void) {return 1; }
57#endif
58
59// // TODO: DUE to the use of HALT in npolygon.cc :(((
60extern "C" {void m2_end(int i){exit(i);}}
61
62// // TODO: DUE to its use in kutil.cc :(((
63char * showOption(){return NULL;}
64
65// // TODO: DUE to its use in feread.cc :(((
66char *iiArithGetCmd(int nPos){return NULL; }
67
68
69#include <coeffs/numbers.h>
70
71#include "structs.h"
72
73
74// HEADERS:
75#include "hutil.h"
76//#include "idrec.h" // moved to Singular
77#include "stairc.h"
78#include "ideals.h"
79#include "syz.h"
80#include "fast_maps.h"
81#include "febase.h"
82#include "walkProc.h"
83#include "walkMain.h"
84#include "walkSupport.h"
85#include "khstd.h"
86/// #include "sparsmat.h" // TODO: install polys/this!
87//+
88
89#include "fglm.h"
90#include "kstd1.h"
91#include "fglmgauss.h"
92#include "fglmvec.h"
93#include "kstdfac.h"
94#include "kmatrix.h"
95#include "GMPrat.h"
96#include "multicnt.h"
97#include "npolygon.h"
98#include "semic.h"
99#include "spectrum.h"
100#include "splist.h"
101#include "multicnt.h"
102#include "eigenval.h"
103#include "units.h"
104#include "ratgring.h"
105#include "shiftgb.h"
106// #include "mmalloc.h" // move to Singular!?
107#include "gfan.h"
108
109#include "kutil.h"
110
111// #include "dbm_sl.h" // TODO: needs si_link// already moved to Singular/!
112
113// #include "CCRing.h" // Too old!
114#include "digitech.h"
115#include "eigenval.h"
116#include "fast_maps.h"
117#include "fast_mult.h"
118#include "febase.h"
119
120#include "fglmgauss.h"
121#include "fglm.h"
122#include "fglmvec.h"
123
124////////#include "F5cData.h"
125#include "f5c.h"
126#include "f5data.h"
127#include "f5gb.h"
128#include "f5lists.h"
129////////#include "F5cLists.h"
130
131#include "gfan.h"
132
133
134#include "GMPrat.h"
135
136// #include "htmlhelp.h" // For Windows //
137#include "hutil.h"
138// #include "Ideal.h" // Too old?
139
140
141#include "ideals.h"
142
143#include "kmatrix.h"
144#include "kstd1.h"
145#include "kstdfac.h"
146#include "khstd.h"
147
148#include "linearAlgebra.h"
149
150
151
152// #include "lplist.h" // Too old!
153#include "multicnt.h"
154#include "npolygon.h"
155// #include "Number.h" // Too old?
156#include "page.h"
157// #include "Poly.h" // Too old?
158// #include "PowerSeries.h" // Too old?
159
160#include "preimage.h"
161
162#include "nc.h"
163
164#include "ratgring.h"
165#include "ringgb.h"
166// #include "run.h" // moved to Singular
167#include "semic.h"
168#include "shiftgb.h"
169// #include "si_gmp.h" // remove it...
170#include "spectrum.h"
171#include "splist.h"
172#include "stairc.h"
173#include "structs.h"
174#include "syz.h"
175// #include "testpoly.h" // Too old?
176
177#include "tgbgauss.h"
178#include "tgb.h"
179
180#include "timer.h"
181
182#include "units.h"
183#include "walkMain.h"
184#include "walkProc.h"
185#include "walkSupport.h"
186
187
188// #include <polys/clapconv.h> // due to factory? :(
189// #include "tgb_internal.h" // :(
190// #include "F4.h" // uses tgb_internal // :(
191// #include "IIntvec.h" // :(
192
193
194
195
196// #include "fglmzero.cc" // looks like <factory/templates/ftmpl_list.h> must be installed!
197// TODO: looks like <coeffs/mpr_complex.h> must be installed!
198
199
200
201#include "polys.h"
202
203void TestGBEngine()
204{
205
206  //  R = MPolynomialRing_polydict(QQ,5,'w,x,y,z,C', order='degrevlex')
207  //  J = (w*w - x*z, w*x - y*z, x*x - w*y, x*y - z*z, y*y - w*z)
208
209  const short w = 1;
210  const short x = 2;
211  const short y = 3;
212  const short z = 4;
213
214  const short N = (z - w + 1);
215
216  char **n=(char**)omalloc(N*sizeof(char*));
217
218
219  n[w-1]=omStrDup("w");
220  n[x-1]=omStrDup("x");
221  n[y-1]=omStrDup("y");
222  n[z-1]=omStrDup("z");
223
224
225  const int D = 3;
226  int *order = (int *) omAlloc0(D* sizeof(int));
227  int *block0 = (int *)omAlloc0(D * sizeof(int));
228  int *block1 = (int *)omAlloc0(D * sizeof(int));
229
230  order[0]  = ringorder_dp;
231  block0[0] = 1;
232  block1[0] = N;
233
234  order[1]  = ringorder_C;
235  block0[1] = 1;
236  block1[1] = N;
237
238  ring R = rDefault(0, N, n, D, order, block0, block1);
239
240//   ring R = rDefault(0, N, n);
241
242  rWrite(R); PrintLn();
243
244#ifdef RDEBUG
245  rDebugPrint(R);
246#endif
247
248  ideal I = idInit(5, 1);
249
250  int gen = 0;
251
252  {
253    // -xz
254    poly p = p_ISet(-1,R);
255
256    p_SetExp(p, x, 1, R);
257    p_SetExp(p, z, 1, R);
258    p_Setm(p, R);
259
260    assume( p_GetExp(p, x, R) == 1 );
261    assume( p_GetExp(p, z, R) == 1 );
262    assume( p_GetExp(p, w, R) == 0 );
263    assume( p_GetExp(p, y, R) == 0 );
264
265    // +w2
266    poly lp = p_ISet(1,R);
267    p_SetExp(lp, w, 2, R);
268    p_Setm(lp, R);
269
270    assume( p_GetExp(lp, w, R) == 2 );
271    assume( p_GetExp(lp, x, R) == 0 );
272    assume( p_GetExp(lp, y, R) == 0 );
273    assume( p_GetExp(lp, z, R) == 0 );
274
275    MATELEM(I, 1, ++gen) = p_Add_q(lp, p, R); // w2 - xz
276  }
277
278  {
279    // -yz
280    poly p = p_ISet(-1,R);
281
282    p_SetExp(p, y, 1, R);
283    p_SetExp(p, z, 1, R);
284    p_Setm(p, R);
285
286    assume( p_GetExp(p, y, R) == 1 );
287    assume( p_GetExp(p, z, R) == 1 );
288    assume( p_GetExp(p, w, R) == 0 );
289    assume( p_GetExp(p, x, R) == 0 );
290
291    // +wx
292    poly lp = p_ISet(1,R);
293    p_SetExp(lp, w, 1, R);
294    p_SetExp(lp, x, 1, R);
295    p_Setm(lp, R);
296
297    assume( p_GetExp(lp, w, R) == 1 );
298    assume( p_GetExp(lp, x, R) == 1 );
299    assume( p_GetExp(lp, y, R) == 0 );
300    assume( p_GetExp(lp, z, R) == 0 );
301
302    MATELEM(I, 1, ++gen) = p_Add_q(lp, p, R); // wx - yz
303  }
304
305
306  {
307    // -wy
308    poly p = p_ISet(-1,R);
309
310    p_SetExp(p, y, 1, R);
311    p_SetExp(p, w, 1, R);
312    p_Setm(p, R);
313
314    assume( p_GetExp(p, y, R) == 1 );
315    assume( p_GetExp(p, w, R) == 1 );
316    assume( p_GetExp(p, z, R) == 0 );
317    assume( p_GetExp(p, x, R) == 0 );
318
319    // +x2
320    poly lp = p_ISet(1,R);
321    p_SetExp(lp, x, 2, R);
322    p_Setm(lp, R);
323
324    assume( p_GetExp(lp, w, R) == 0 );
325    assume( p_GetExp(lp, x, R) == 2 );
326    assume( p_GetExp(lp, y, R) == 0 );
327    assume( p_GetExp(lp, z, R) == 0 );
328
329    MATELEM(I, 1, ++gen) = p_Add_q(lp, p, R); // x2 - wy
330  }
331
332
333  {
334    // -z2
335    poly p = p_ISet(-1,R);
336
337    p_SetExp(p, z, 2, R);
338    p_Setm(p, R);
339
340    assume( p_GetExp(p, y, R) == 0 );
341    assume( p_GetExp(p, w, R) == 0 );
342    assume( p_GetExp(p, z, R) == 2 );
343    assume( p_GetExp(p, x, R) == 0 );
344
345    // +xy
346    poly lp = p_ISet(1,R);
347    p_SetExp(lp, x, 1, R);
348    p_SetExp(lp, y, 1, R);
349    p_Setm(lp, R);
350
351    assume( p_GetExp(lp, w, R) == 0 );
352    assume( p_GetExp(lp, x, R) == 1 );
353    assume( p_GetExp(lp, y, R) == 1 );
354    assume( p_GetExp(lp, z, R) == 0 );
355
356    MATELEM(I, 1, ++gen) = p_Add_q(lp, p, R); // xy - z2
357  }
358
359
360  {
361    // -wz
362    poly p = p_ISet(-1,R);
363
364    p_SetExp(p, w, 1, R);
365    p_SetExp(p, z, 1, R);
366    p_Setm(p, R);
367
368    assume( p_GetExp(p, y, R) == 0 );
369    assume( p_GetExp(p, w, R) == 1 );
370    assume( p_GetExp(p, z, R) == 1 );
371    assume( p_GetExp(p, x, R) == 0 );
372
373    // +y2
374    poly lp = p_ISet(1,R);
375    p_SetExp(lp, y, 2, R);
376    p_Setm(lp, R);
377
378    assume( p_GetExp(lp, w, R) == 0 );
379    assume( p_GetExp(lp, x, R) == 0 );
380    assume( p_GetExp(lp, y, R) == 2 );
381    assume( p_GetExp(lp, z, R) == 0 );
382
383    MATELEM(I, 1, ++gen) = p_Add_q(lp, p, R); // y2 - wz
384  }
385#ifdef PDEBUG
386  PrintS("I: ");
387  idShow(I, R, R, 0);
388#endif
389
390
391//  ideal kStd(ideal F, ideal Q, tHomog h, intvec ** mw,intvec *hilb=NULL,
392//             int syzComp=0,int newIdeal=0, intvec *vw=NULL);
393  // make R the default ring:
394  rChangeCurrRing(R);
395
396  {
397    ideal G = kStd(I, currQuotient, testHomog, NULL);
398
399#ifdef PDEBUG
400    PrintS("GB: ");
401    idShow(G, R, R, 0);
402#endif
403
404    idDelete( &G, R);
405  }
406
407  {
408    intvec *weights = NULL;
409    ideal SYZ = idSyzygies(I, testHomog, &weights);
410
411#ifdef PDEBUG
412    PrintS("SYZ: ");
413    idShow(SYZ, R, R, 0);
414#endif
415
416    idDelete(&SYZ, R);
417    if (weights!=NULL) { PrintS("weights: "); weights->show(); delete weights; }
418  }
419
420
421  {
422    PrintS("\n**********************************\n");
423    PrintS("lres: \n");
424    int dummy;
425    syStrategy r = syLaScala3(I,&dummy);
426
427    intvec *b = syBettiOfComputation(r, FALSE);
428    PrintS("non-min. betti: \n");    b->show();    PrintLn();
429    delete b;
430
431    Print("length: %d\n", sySize(r));
432
433    syPrint(r, "R");
434
435    r =  syMinimize(r); // syzstr->references ++ ==> memory leak :(((
436
437    b = syBettiOfComputation(r, TRUE);
438    PrintS("min. betti: \n");    b->show();    PrintLn();
439    delete b;
440
441    Print("length: %d\n", sySize(r));
442
443    syPrint(r, "R");
444
445    syKillComputation(r, R);
446  }
447
448  {
449    PrintS("\n**********************************\n");
450    PrintS("sres: \n");
451    const int maxl = rVar(R)-1; // +2*(1);
452
453    syStrategy r = sySchreyer(I, rVar(R));
454
455    intvec *b = syBettiOfComputation(r, FALSE);
456    PrintS("non-min. betti: \n");    b->show();    PrintLn();
457    delete b;
458
459    Print("length: %d\n", sySize(r));
460
461    syPrint(r, "R");
462
463    r =  syMinimize(r); // syzstr->references ++ ==> memory leak :(((
464
465    b = syBettiOfComputation(r, TRUE);
466    PrintS("min. betti: \n");    b->show();    PrintLn();
467    delete b;
468
469    Print("length: %d\n", sySize(r));
470
471    syPrint(r, "R");
472
473    syKillComputation(r, R);
474  }
475
476
477
478  {
479    PrintS("\n**********************************\n");
480    PrintS("nres: \n");
481    intvec *weights=NULL;
482//    const int maxl = rVar(R)-1 + 2*(1);
483    syStrategy r = syResolution(I, rVar(R)-1, weights, FALSE/*iiOp==MRES_CMD*/);
484
485    intvec *b = syBettiOfComputation(r, FALSE);
486    PrintS("non-min. betti: \n");    b->show();    PrintLn();
487    delete b;
488
489    Print("length: %d\n", sySize(r));
490
491    syPrint(r, "R");
492
493    r =  syMinimize(r); // syzstr->references ++ ==> memory leak :(((
494
495    b = syBettiOfComputation(r, TRUE);
496    PrintS("min. betti: \n");    b->show();    PrintLn();
497    delete b;
498
499    Print("length: %d\n", sySize(r));
500
501    syPrint(r, "R");
502
503    syKillComputation(r, R);
504  }
505
506
507  {
508    PrintS("\n**********************************\n");
509    PrintS("mres: \n");
510    intvec *weights=NULL;
511//    const int maxl = rVar(R)-1 + 2*(1);
512    syStrategy r = syResolution(I, rVar(R)+1, weights, TRUE/*iiOp==MRES_CMD*/);
513
514    intvec *b = syBettiOfComputation(r, FALSE);
515    PrintS("non-min. betti: \n");    b->show();    PrintLn();
516    delete b;
517
518    Print("length: %d\n", sySize(r));
519
520    syPrint(r, "R");
521
522    r =  syMinimize(r); // syzstr->references ++ ==> memory leak :(((
523
524    b = syBettiOfComputation(r, TRUE);
525    PrintS("min. betti: \n");    b->show();    PrintLn();
526    delete b;
527
528    Print("length: %d\n", sySize(r));
529
530    syPrint(r, "R");
531
532    syKillComputation(r, R);
533  }
534
535
536
537
538  idDelete( &I, R);
539  rDelete(R); // should cleanup every belonging polynomial, right!?
540
541}
542
543
544
545void TestSimpleRingArithmetcs()
546{
547  // Libpolys tests:
548
549  // construct the ring Z/32003[x,y,z]
550  // the variable names
551  char **n=(char**)omalloc(3*sizeof(char*));
552  n[0]=omStrDup("x");
553  n[1]=omStrDup("y");
554  n[2]=omStrDup("z2");
555
556  ring R = rDefault(32003,3,n); //  ring R = rDefault(0,3,n);
557
558  rWrite(R); PrintLn();
559
560#ifdef RDEBUG
561  rDebugPrint(R);
562#endif
563
564
565  poly p = p_ISet(1,R); p_SetExp(p,1,1, R); p_Setm(p, R);
566
567  assume( p_GetExp(p,1, R) == 1 );
568
569  poly pp = pp_Mult_qq( p, p, R);
570
571  Print("p: "); p_Write0(p, R); Print(", deg(p): %d", p_Totaldegree(p, R)); assume( 1 == p_Totaldegree(p, R) );
572
573  Print("; p*p : "); p_Write0(pp, R); Print("deg(pp): %d\n", p_Totaldegree(pp, R)); assume( 2 == p_Totaldegree(pp, R) );
574
575
576  p_Delete(&p, R);
577
578  assume( p_GetExp(pp,1, R) == 2 );
579
580  p_Delete(&pp, R);
581
582
583//  rDelete(R);
584
585  // make R the default ring:
586  rChangeCurrRing(R);
587
588  // create the polynomial 1
589  poly p1=pISet(1);
590
591  // create tthe polynomial 2*x^3*z^2
592  poly p2=p_ISet(2,R);
593  pSetExp(p2,1,3);
594  pSetExp(p2,3,2);
595  pSetm(p2);
596
597  // print p1 + p2
598  Print("p1: "); pWrite0(p1);
599  Print(" + p2: "); pWrite0(p2);
600  Print("  ---- >>>> ");
601
602  // compute p1+p2
603  p1=p_Add_q(p1,p2,R); p2=NULL;
604  pWrite(p1);
605
606  // clean up:
607//  pDelete(&p1);
608
609  rDelete(R); // should cleanup every belonging polynomial, right!?
610}
611
612
613int main( int, char *argv[] )
614{
615  feInitResources(argv[0]);
616
617  StringSetS("ressources in use (as reported by feStringAppendResources(0):\n");
618  feStringAppendResources(0);
619  PrintS(StringAppendS("\n"));
620
621  TestGBEngine();
622  TestSimpleRingArithmetcs();
623
624  return 0;
625}
Note: See TracBrowser for help on using the repository browser.