Annotation of OpenXM/src/kan96xx/Doc/gfan.sm1, Revision 1.10
1.10 ! takayama 1: % $OpenXM: OpenXM/src/kan96xx/Doc/gfan.sm1,v 1.9 2005/06/30 08:39:39 takayama Exp $
1.1 takayama 2: % cp cone.sm1 $OpenXM_HOME/src/kan96xx/Doc/gfan.sm1
1.10 ! takayama 3: % $Id: cone.sm1,v 1.74 2005/07/07 01:32:07 taka Exp $
1.1 takayama 4: % iso-2022-jp
1.9 takayama 5: %%Ref: @s/2004/08/21-note.pdf
1.1 takayama 6:
1.6 takayama 7: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
8: %% Two examples are given below to get a global Grobner fan and
9: %% a local Grobner fan ; cone.sample and cone.sample2
10: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
11: %%% Global Grobner Fan
12: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
13: %% How to input data? An example. (cf. test13.sm1)
14: %% Modify the following or copy the /cone.sample { ... } def
15: %% to your own file,
1.9 takayama 16: %% edit it, and execute it by " cone.sample ; "
1.6 takayama 17: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
18: /cone.sample {
19: cone.load.cohom
20: % write a comment about the problem. "nl" means new line.
21: /cone.comment [
22: (Toric ideal for 1-simplex x 2-simplex, in k[x]) nl
23: ] cat def
24:
25: % List of variables
26: % If cone.type=1, then (H) should be added.
27: /cone.vlist [(x11) (x12) (x13) (x21) (x22) (x23)
28: (Dx11) (Dx12) (Dx13) (Dx21) (Dx22) (Dx23) (h)] def
29:
30: % List of variables in the form for define_ring.
31: /cone.vv (x11,x12,x13,x21,x22,x23) def
32:
33: % If cone.type=0, then x,Dx,
34: % If cone.type=1, then x,Dx,h,H (Doubly homogenized)
35: % If cone.type=2, then x,Dx,h
36: /cone.type 2 def
37:
38: % Set how to parametrize the weight space.
39: % In the example below, 6 means the number of variables x11,x12,x13,x21,x22,x33
40: % p q parametrizeSmallFan (p >= q) : Enumerate Grobner cones in the Small
41: % Grobner fan.
42: % The weights for the last p-q variables
43: % are 0.
44: % Example. 6 2 parametrizeSmallFan weights for x12,x21,x22,x23 are 0.
45: %
46: % p q parametrizeTotalFan (p = q = number of variables in cone.vv)
47: % p > q has not yet been implemented.
48: %
49: /cone.parametrizeWeightSpace {
50: 6 6 parametrizeSmallFan
51: } def
52:
53: % If you want to enumerate Grobner cones in local order (i.e., x^e <= 0),
54: % then cone.local = 1 else cone.local = 0.
55: /cone.local 0 def
56:
57: % Initial value of the weight in the weight space of which dimension is
58: % cone.m
59: % If it is null, then a random weight is used.
60: /cone.w_start
61: null
62: def
63:
64: % If cone.h0=1, then the weight for h is 0.
65: % It is usally set to 1.
66: /cone.h0 1 def
67:
68: % Set input polynomials which generate the ideal.
69: % Input must be homogenized.
70: % (see also data/test14.sm1 for double homogenization.)
71: /cone.input
72: [
73: (x11 x22 - x12 x21)
74: (x12 x23 - x13 x22)
75: (x11 x23 - x13 x21)
76: ]
77: def
78:
1.10 ! takayama 79: /cone.DhH 0 def
1.6 takayama 80: % Set a function to compute Grobner basis.
81: % cone.gb_Dh : For computing in Homogenized Weyl algebra h[1,1](D).
82: % cone.gb_DhH : For computing in doubly homogenized Weyl algebra.
83: % ( Computation in ^O and h[0,1](^D) need this
84: % as the first step. /cone.local 1 def )
85: /cone.gb {
86: cone.gb_Dh
87: } def
88:
89:
90: cone.comment message
91: (cone.input = ) message
92: cone.input message
93: %%%% Step 1. Enumerating the Grobner Cones in a global ring.
94: %%%% The result is stored in cone.fan
95: getGrobnerFan
96:
97: %%%% If you want to print the output, then uncomment.
98: printGrobnerFan
99:
100: %%%% If you want to save the data to the file sm1out.txt, then uncomment.
1.9 takayama 101: % /cone.withGblist 1 def saveGrobnerFan /ff set ff output
1.6 takayama 102:
103: %%%% Step 2. Dehomogenize the Grobner Cones
104: %%%% by the equivalence relation in a local ring (uncomment).
105: % dhCones_h
106:
107: %%%% Generate the final data dhcone2.fan (a list of local Grobner cones.)
108: % dhcone.rtable
109:
110: %%%% Output dhcone2.fan with explanations
111: % dhcone.printGrobnerFan
112:
113: } def
114: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
115: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
116: %% End of " How to input data? An example. "
117: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
118: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
119:
120:
121:
122: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
123: %%% Local Grobner Fan
124: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
125: %% How to input data? The example 2 (cf. test14.sm1).
126: %% Modify the following or copy the /cone.sample2 { ... } def
127: %% to your own file,
128: %% edit it, and execute if by " cone.sample2 ; "
129: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
130: /cone.sample2 {
131: cone.load.cohom
132: % write a comment about the problem. "nl" means new line.
133: /cone.comment [
134: (BS for y and y-(x-1)^2, t1, t2 space, in doubly homogenized Weyl algebra.) nl
135: (The Grobner cones are dehomogenized to get local Grobner fan.) nl
136: ] cat def
137:
138: % List of variables
139: % If cone.type=1, then (H) should be added.
140: /cone.vlist [(t1) (t2) (x) (y) (Dt1) (Dt2) (Dx) (Dy) (h) (H)] def
141:
142: % List of variables in the form for define_ring.
143: /cone.vv (t1,t2,x,y) def
144:
145: % If cone.type=0, then x,Dx,
146: % If cone.type=1, then x,Dx,h,H (Doubly homogenized)
147: % If cone.type=2, then x,Dx,h
148: /cone.type 1 def
149:
150: % Set how to parametrize the weight space.
151: % In the example below, 6 means the number of variables x11,x12,x13,x21,x22,x33
152: % p q parametrizeSmallFan (p >= q) : Enumerate Grobner cones in the Small
153: % Grobner fan.
154: % The weights for the last p-q variables
155: % are 0.
156: % Example. 6 2 parametrizeSmallFan weights for x12,x21,x22,x23 are 0.
157: %
158: % p q parametrizeTotalFan (p = q = number of variables in cone.vv)
159: % p > q has not yet been implemented.
160: %
161: /cone.parametrizeWeightSpace {
162: 4 2 parametrizeSmallFan
163: } def
164:
165: % If you want to enumerate Grobner cones in local order (i.e., x^e <= 0),
166: % then cone.local = 1 else cone.local = 0.
167: /cone.local 1 def
168:
169: % Initial value of the weight in the weight space of which dimension is
170: % cone.m
171: % If it is null, then a random weight is used.
172: /cone.w_start
173: null
174: def
175:
176: % If cone.h0=1, then the weight for h is 0.
177: % It is usally set to 1.
178: /cone.h0 1 def
179:
180: % Set input polynomials which generate the ideal.
181: % Input must be homogenized.
182: % (see also data/test14.sm1 for double homogenization.)
183: /cone.input
184: [
185: (t1-y) (t2 - (y-(x-1)^2))
186: ((-2 x + 2)*Dt2+Dx)
187: (Dt1+Dt2+Dy)
188: ]
189: def
190: % homogenize
191: [cone.vv ring_of_differential_operators
192: [[(t1) -1 (t2) -1 (Dt1) 1 (Dt2) 1]] ecart.weight_vector
193: 0] define_ring
194: dh.begin
195: cone.input { . homogenize toString } map /cone.input set
196: dh.end
197:
1.10 ! takayama 198: /cone.DhH 1 def
1.6 takayama 199: % Set a function to compute Grobner basis.
200: % cone.gb_Dh : For computing in Homogenized Weyl algebra h[1,1](D).
201: % cone.gb_DhH : For computing in doubly homogenized Weyl algebra.
202: % ( Computation in ^O and h[0,1](^D) need this
203: % as the first step. /cone.local 1 def )
204: /cone.gb {
205: cone.gb_DhH
206: } def
207:
208: cone.comment message
209: (cone.input = ) message
210: cone.input message
211: %%%% Step 1. Enumerating the Grobner Cones in a global ring.
212: %%%% The result is stored in cone.fan
213: getGrobnerFan
214:
215: %%%% If you want to print the output, then uncomment.
216: printGrobnerFan
217:
218: %%%% If you want to save the data to the file sm1out.txt, then uncomment.
1.9 takayama 219: % /cone.withGblist 1 def saveGrobnerFan /ff set ff output
1.6 takayama 220:
221: %%%% Step 2. Dehomogenize the Grobner Cones
222: %%%% by the equivalence relation in a local ring (uncomment).
223: dhCones_h
224:
225: %%%% Generate the final data dhcone2.fan (a list of local Grobner cones.)
226: dhcone.rtable
227:
228: %%%% Output dhcone2.fan with explanations
229: dhcone.printGrobnerFan
230:
231: } def
232: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
233: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
234: %% End of " How to input data? The example 2. "
235: %% Do not touch below.
236: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
237: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
238:
239:
240: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
241: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
242: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
243: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
244:
245: [(parse) (cgi.sm1) pushfile] extension
246:
247: % If you use local polymake, then comment out.
248: % If you use the cgi/polymake on the net, then uncomment out.
1.8 takayama 249: %/doPolymake {doPolymake.OoHG} def (Using doPolymake.OoHG ) message
250: %/polymake.start {polymake.start.OoHG} def (Using polymake.start.OoHG ) message
251: %% Choose it automatically.
252: [(which) (polymake)] oxshell tag 0 eq {
253: (Polymake is not installed in this system.) message
254: /doPolymake {doPolymake.OoHG} def
255: (Using doPolymake.OoHG ) message
256: /polymake.start {polymake.start.OoHG} def
257: (Using polymake.start.OoHG ) message
258: } { (Local polymake will be used.) message } ifelse
1.6 takayama 259:
1.1 takayama 260: /cone.debug 1 def
261:
262: /ox.k0.loaded boundp {
263: } {
264: [(parse) (ox.sm1) pushfile] extension
265: } ifelse
266:
1.6 takayama 267: /cone.load.cohom {
268: /cone.loaded boundp { }
269: {
270: [(parse) (cohom.sm1) pushfile] extension
1.8 takayama 271: % [(parse) (cone.sm1) pushfile] extension % BUG? cone.sm1 overrides a global
272: % in cohom.sm1?
1.6 takayama 273: [(parse) (dhecart.sm1) pushfile] extension
274: /cone.loaded 1 def
1.8 takayama 275: oxNoX
276: polymake.start ( ) message
1.6 takayama 277: } ifelse
278: } def
279:
280: %% Usages: cone.gb_DhH. h H (double homogenized) $BMQ$N(B GB.
281: %% dhecart.sm1 $B$r(B load $B$7$F$"$k$3$H(B. $BF~NO$OF1<!$G$J$$$H$$$1$J$$(B.
282: %% [cone.vv ring_of_differential_operators
283: %% [[(t1) -1 (t2) -1 (Dt1) 1 (Dt2) 1]] ecart.weight_vector
284: %% 0] define_ring
285: %% dh.begin homogenize dh.end $B$J$I$NJ}K!$GF1<!2=$G$-$k(B.
286: /cone.gb_DhH {
287: /arg2 set /arg1 set
288: [/ff /ww] pushVariables
289: [
290: /ff arg1 def
291: /ww arg2 def
292: /dh.gb.verbose 1 def
293: /dh.autoHomogenize 0 def
294: [(AutoReduce) 1] system_variable
295: [ff { toString } map cone.vv
1.9 takayama 296: [ww cone.vv generateD1_1]] ff getAttributeList setAttributeList
297: dh.gb 0 get /arg1 set
1.6 takayama 298: ] pop
1.9 takayama 299: popVariables
1.6 takayama 300: arg1
301: } def
302:
1.3 takayama 303: %
304: % cone.fan, cone.gblist $B$K(B fan $B$N%G!<%?$,$O$$$k(B.
305: %
1.6 takayama 306: %%%%<<<< $B=i4|%G!<%?$N@_DjNc(B. $BF|K\8lHG(B data/test13 $B$h$j(B. <<<<<<<<<<<<<<
307: /cone.sample.test13.ja {
1.2 takayama 308: /cone.loaded boundp { }
309: {
310: [(parse) (cohom.sm1) pushfile] extension
311: [(parse) (cone.sm1) pushfile] extension
312: /cone.loaded 1 def
313: } ifelse
314: /cone.comment [
315: (Toric ideal for 1-simplex x 2-simplex, in k[x]) nl
316: ] cat def
317: %------------------Globals----------------------------------------
318: % Global: cone.type
319: % $B$I$N(B exponents $B$r<h$j=P$9$N$+;XDj$9$k(B.
320: % cf. exponents, gbext h $B$d(B H $B$b8+$k$+(B?
321: % 0 : x,y,Dx,Dy
322: % 1 : x,y,Dx,Dy,h,H
323: % 2 : x,y,Dx,Dy,h
324: /cone.type 2 def
325:
326: % Global: cone.local
327: % cone.local: Local $B$+(B? 1 $B$J$i(B local
328: /cone.local 0 def
329:
330:
331: % Global: cone.h0
332: % cone.h0: 1 $B$J$i(B h $B$N(B weight 0 $B$G$N(B Grobner fan $B$r7W;;$9$k(B.
333: /cone.h0 1 def
334:
335: % --------------- $BF~NO%G!<%?MQBg0hJQ?t$N@_Dj(B --------------------------
336: %
337: % cone.input : $BF~NOB?9`<07O(B
338: /cone.input
339: [
340: (x11 x22 - x12 x21) (x12 x23 - x13 x22)
341: (x11 x23 - x13 x21)
342: ]
343: def
344:
345: % cone.vlist : $BA4JQ?t$N%j%9%H(B
346: /cone.vlist [(x11) (x12) (x13) (x21) (x22) (x23)
347: (Dx11) (Dx12) (Dx13) (Dx21) (Dx22) (Dx23) (h)] def
348:
349: % cone.vv : define_ring $B7A<0$NJQ?t%j%9%H(B.
350: /cone.vv (x11,x12,x13,x21,x22,x23) def
351:
352: % cone.parametrizeWeightSpace : weight $B6u4V$r(B parametrize $B$9$k4X?t(B.
353: % $BBg0hJQ?t(B cone.W , cone.Wpos $B$b$-$^$k(B.
354: /cone.parametrizeWeightSpace {
355: 6 6 parametrizeSmallFan
356: } def
357:
358: % cone.w_start : weight$B6u4V$K$*$1$k(B weight $B$N=i4|CM(B.
359: % $B$3$NCM$G(B max dim cone $B$,F@$i$l$J$$$H(B random weight $B$K$h$k(B $B%5!<%A$,;O$^$k(B.
360: % random $B$K$d$k$H$-$O(B null $B$K$7$F$*$/(B.
361: /cone.w_start
362: [9 8 5 4 5 6]
363: def
364:
365: % cone.gb : gb $B$r7W;;$9$k4X?t(B.
366: /cone.gb {
367: cone.gb_Dh
368: } def
369:
370:
371:
372: ( ) message
373: cone.comment message
374: (cone.input = ) messagen cone.input message
375: (Type in getGrobnerFan) message
376: (Do clearGlobals if necessary) message
377: (printGrobnerFan ; saveGrobnerFan /ff set ff output ) message
378:
379: } def
380: %%%%%%>>>>> $B=i4|%G!<%?$N@_DjNc$*$o$j(B >>>>>>>>>>>>>>>>>>>>>>
381:
1.1 takayama 382: % Global: cone.type
383: % $B$I$N(B exponents $B$r<h$j=P$9$N$+;XDj$9$k(B.
384: % cf. exponents, gbext h $B$d(B H $B$b8+$k$+(B?
385: % 0 : x,y,Dx,Dy
386: % 1 : x,y,Dx,Dy,h,H
387: % 2 : x,y,Dx,Dy,h
388: /cone.type 2 def
389:
390: % Global: cone.local
391: % cone.local: Local $B$+(B? 1 $B$J$i(B local
392: /cone.local 1 def
393:
394: % Global: cone.h0
395: % cone.h0: 1 $B$J$i(B h $B$N(B weight 0 $B$G$N(B Grobner fan $B$r7W;;$9$k(B.
396: /cone.h0 1 def
397:
398: % Global: cone.n (number of variables in GB)
399: % cone.m (freedom of the weight space. cf. cone.W)
400: % cone.d (pointed cones lies in this space. cf. cone.Lp)
401: % These are set during getting the cone.startingCone
402:
1.10 ! takayama 403: %<
! 404: % global
! 405: %cone.ckmFlip. Collar-Kalkbrener-Mall $B$N(B flip $B%"%k%4%j%:%`$r;H$o$J$$(B 0. $B;H$&(B 1.
! 406: % Default $B$O(B 0.
! 407: %>
! 408: /cone.ckmFlip 0 def
! 409:
! 410: %<
! 411: % global
! 412: % cone.DhH dx x = x dx + h H $B$J$i(B 1. dx x = x dx + h^2 $B$J$i(B 0. Default 0.
! 413: %>
! 414: /cone.DhH 0 def
! 415:
! 416: % Default $B$N(B cone.gb $B$NDj5A(B. $B3F%W%m%0%i%`$G:FEYDj5A$7$F$b$h$$(B.
! 417: /cone.gb {
! 418: cone.DhH {
! 419: cone.gb_DhH
! 420: } {
! 421: cone.gb_Dh
! 422: } ifelse
! 423: } def
1.1 takayama 424:
425: %<
426: % Usage: wv g coneEq1
427: % in(f) $B$,(B monomial $B@lMQ(B. in_w(f) = LT(f) $B$H$J$k(B weight w $B$NK~$?$9(B
428: % $BITEy<0@)Ls$r5a$a$k(B.
429: %>
430: /coneEq1 {
431: /arg1 set
432: [/g /eqs /gsize /i /j /n /f /exps /m % Do not use "eq" as a variable
433: /expsTop
434: ] pushVariables
435: [
436: /g arg1 def % Reduced Grobner basis
437: /eqs [ ] def % $BITEy<07O$N78?t(B
438: /gsize g length def
439: 0 1 gsize 1 sub {
440: /i set
441: g i get /f set % f $B$O(B i $BHVL\$N(B reduced Grobner basis $B$N85(B
442: [(exponents) f cone.type] gbext /exps set % exps $B$O(B f $B$N(B exponent vector
443: exps length /m set
444: m 1 eq not {
445: /expsTop exps 0 get def % expsTop $B$O(B f $B$N@hF,$N(B exponent vector.
446: 1 1 exps length 1 sub {
447: /j set
448: eqs [expsTop exps j get sub] join /eqs set
449: % exps[0]-exps[j] $B$r(B eqs $B$X3JG<$7$F$$$/$@$1(B.
450: % Cone $B$N(B closure $B$r$@$9$N$G(B >= $B$G(B OK.
451: } for
452: } { } ifelse
453: } for
454: /arg1 eqs def
455: ] pop
456: popVariables
457: arg1
458: } def
459:
460: %<
461: % Usage: ww g coneEq
462: % ww $B$O(B [v1 w1 v2 w2 ... ] $B7A<0(B. (v-w $B7A<0(B) w1, w2 $B$O(B univNumber $B$G$b$$$$(B.
463: % g $B$O(B reduced Grobner basis
464: % in(f) $B$,(B monomial $B$G$J$$>l9g$b07$&(B.
465: % in_w(f) = in_ww(f) $B$H$J$k(B weight w $B$NK~$?$9(B
466: % $BITEy<0@)Ls$r5a$a$k(B.
467: % ord_w, init (weightv) $B$rMQ$$$k(B.
468: %>
469: /coneEq {
470: /arg2 set
471: /arg1 set
472: [/g /eqs /gsize /i /j /n /f /exps /m
473: /expsTop /ww /ww2 /iterms
474: ] pushVariables
475: [
476: /g arg2 def % Reduced Grobner basis
477: /ww arg1 def % weight vector. v-w $B7A<0(B
478: ww to_int32 /ww set % univNum $B$,$"$l$P(B int32 $B$KD>$7$F$*$/(B.
479: /ww2 ww weightv def % v-w $B7A<0$r(B $B?t;z$N%Y%/%H%k$K(B. (init $BMQ(B)
480:
1.3 takayama 481: /eqs null def % $BITEy<07O$N78?t(B
1.1 takayama 482: /gsize g length def
483: 0 1 gsize 1 sub {
484: /i set
485: g i get /f set % f $B$O(B i $BHVL\$N(B reduced Grobner basis $B$N85(B
486: [(exponents) f cone.type] gbext /exps set % exps $B$O(B f $B$N(B exponent vector
487: exps length /m set
488: m 1 eq not {
489: /expsTop exps 0 get def % expsTop $B$O(B f $B$N@hF,$N(B exponent vector.
490: /iterms f ww2 init length def % f $B$N(B initial term $B$N9`$N?t(B.
491: % in_ww(f) > f_j $B$H$J$k9`$N=hM}(B.
492: iterms 1 exps length 1 sub {
493: /j set
1.3 takayama 494: expsTop exps j get sub eqs cons /eqs set
1.1 takayama 495: % exps[0]-exps[j] $B$r(B eqs $B$X3JG<$7$F$$$/(B.
496: } for
497: % in_ww(f) = f_j $B$H$J$k9`$N=hM}(B.
498: [(exponents) f ww2 init cone.type] gbext /exps set % exps $B$O(B in(f)
499: 1 1 iterms 1 sub {
500: /j set
1.3 takayama 501: exps j get expsTop sub eqs cons /eqs set
502: expsTop exps j get sub eqs cons /eqs set
1.1 takayama 503: % exps[j]-exps[0], exps[0]-exps[j] $B$r3JG<(B.
504: % $B7k2LE*$K(B (exps[j]-exps[0]).w = 0 $B$H$J$k(B.
505: } for
506: } { } ifelse
507: } for
1.3 takayama 508: eqs listToArray reverse /eqs set
1.1 takayama 509: /arg1 eqs def
510: ] pop
511: popVariables
512: arg1
513: } def
514:
515: %<
516: % Usage: wv g coneEq genPo
517: % polymake $B7A<0$N(B INEQUALITIES $B$r@8@.$9$k(B. coneEq -> genPo $B$HMxMQ(B
518: %>
519: /genPo {
520: /arg1 set
521: [/outConeEq /rr /nn /ii /mm /jj /ee] pushVariables
522: [
523: /outConeEq arg1 def
524: /rr [(INEQUALITIES) nl] cat def % $BJ8;zNs(B rr $B$KB-$7$F$$$/(B.
525: outConeEq length /nn set
526: 0 1 nn 1 sub {
527: /ii set
528: outConeEq ii get /ee set
529: [ rr
530: (0 ) % $BHs$;$$$8MQ$N(B 0 $B$r2C$($k(B.
531: 0 1 ee length 1 sub {
532: /jj set
533: ee jj get toString ( )
534: } for
535: nl
536: ] cat /rr set
537: } for
538: /arg1 rr def
539: ] pop
540: popVariables
541: arg1
542: } def
543:
544: %<
545: % Usage: wv g coneEq genPo2
546: % doPolyamke $B7A<0$N(B INEQUALITIES $B$r@8@.$9$k(B. coneEq -> genPo2 $B$HMxMQ(B
547: % tfb $B7A<0J8;zNs(B.
548: %>
549: /genPo2 {
550: /arg1 set
551: [/outConeEq /rr /nn /ii /mm /jj /ee] pushVariables
552: [
553: /outConeEq arg1 def
554: /rr $polymake.data(polymake.INEQUALITIES([$ def
555: % $BJ8;zNs(B rr $B$KB-$7$F$$$/(B.
556: outConeEq length /nn set
557: 0 1 nn 1 sub {
558: /ii set
559: outConeEq ii get /ee set
560: [ rr
561: ([0,) % $BHs$;$$$8MQ$N(B 0 $B$r2C$($k(B.
562: 0 1 ee length 1 sub {
563: /jj set
564: ee jj get toString
565: jj ee length 1 sub eq { } { (,) } ifelse
566: } for
567: (])
568: ii nn 1 sub eq { } { (,) } ifelse
569: ] cat /rr set
570: } for
571: [rr $]))$ ] cat /rr set
572: /arg1 rr def
573: ] pop
574: popVariables
575: arg1
576: } def
577:
578: /test1 {
579: [(x,y) ring_of_differential_operators 0] define_ring
580: [ (x + y + Dx + Dy).
581: (x ^2 Dx^2 + y^2 Dy^2).
582: (x).
583: ] /gg set
584: gg coneEq1 /ggc set
585: gg message
586: ggc pmat
587:
588: ggc genPo message
589: } def
590:
591: /test2 {
592: [(parse) (dhecart.sm1) pushfile] extension
593: dh.test.p1 /ff set
594: ff 0 get coneEq1 /ggc set
595: ggc message
596: ggc genPo /ss set
597: ss message
598: (Data is in ss) message
599: } def
600:
601:
602: /test3 {
603: % [(parse) (cohom.sm1) pushfile] extension
604: /ww [(Dx) 1 (Dy) 1] def
605: [(x,y) ring_of_differential_operators
606: [ww] weight_vector
607: 0] define_ring
608: [ (x Dx + y Dy -1).
609: (y^2 Dy^2 + 2 + y Dy ).
610: ] /gg set
611: gg {homogenize} map /gg set
612: [gg] groebner 0 get /gg set
613: ww message
614: ww gg coneEq /ggc set
615: gg message
616: ggc pmat
617:
618: ggc genPo message
619: } def
620:
621: %<
622: % Usage: test3b
623: % Grobner cone $B$r7hDj$7$F(B, polymake $BMQ$N%G!<%?$r@8@.$9$k%F%9%H(B.
624: % weight (0,0,1,1) $B$@$H(B max dim cone $B$G$J$$(B.
625: %>
626: /test3b {
627: % [(parse) (cohom.sm1) pushfile] extension
628: /ww [(Dx) 1 (Dy) 2] def
629: [(x,y) ring_of_differential_operators
630: [ww] weight_vector
631: 0] define_ring
632: [ (x Dx + y Dy -1).
633: (y^2 Dy^2 + 2 + y Dy ).
634: ] /gg set
635: gg {homogenize} map /gg set
636: [gg] groebner 0 get /gg set
637: ww message
638: ww gg coneEq /ggc set
639: gg message
640: ggc pmat
641:
642: % ggc genPo /ggs set % INEQ $B$rJ8;zNs7A<0$G(B
643: % ggs message
644: % ggs output
645: % (mv sm1out.txt test3b.poly) system
646: % (Type in polymake-pear.sh test3b.poly FACETS) message
647:
648: ggc genPo2 /ggs set % INEQ $B$rJ8;zNs7A<0(B for doPolymake
649: ggs message
650:
651: } def
652:
653: % commit (dr.sm1): lcm, denominator, ngcd, to_univNum, numerator, reduce
654: % 8/22, changelog-ja $B$^$@(B.
655: % to do : nnormalize_vec, sort_vec --> shell $B$G(B OK.
656: % 8/27, getNode
657:
658: /test4 {
659: $polymake.data(polymake.INEQUALITIES([[0,1,0,0],[0,0,1,0]]))$ /ff set
660: [(FACETS) ff] doPolymake /rr set
661:
662: rr 1 get /rr1 set
663: rr1 getLinearitySubspace pmat
664:
665: } def
666:
667: %<
668: % Usage: vv ineq isInLinearSpace
669: % vv $B$,(B ineq[i] > 0 $B$GDj5A$5$l$kH>6u4V$N$I$l$+$K$O$$$C$F$$$k$J$i(B 0
670: % vv $B$,(B $BA4$F$N(B i $B$K$D$$$F(B ineq[i] = 0 $B$K$O$$$C$F$$$?$i(B 1.
671: %>
672: /isInLinearSpace {
673: /arg2 set
674: /arg1 set
675: [/vv /ineq /ii /rr] pushVariables
676: [
677: /vv arg1 def
678: /ineq arg2 def
679: /rr 1 def
680: {
681: 0 1 ineq length 1 sub {
682: /ii set
683: % vv . ineq[ii] != 0 $B$J$i(B vv $B$O(B linearity space $B$N85$G$J$$(B.
684: vv ineq ii get mul to_univNum isZero {
685: } { /rr 0 def exit} ifelse
686: } for
687: exit
688: } loop
689: /arg1 rr def
690: ] pop
691: popVariables
692: arg1
693: } def
694:
695: %<
696: % Usages: doPolymakeObj getLinearitySubspace
697: % INEQUALITIES $B$H(B VERTICES $B$+$i(B maximal linearity subspace
698: % $B$N@8@.%Y%/%H%k$r5a$a$k(B.
699: % $BNc(B: VERTICES [[0,1,0,0],[0,0,1,0],[0,0,0,-1],[0,0,0,1]]]
700: % $BNc(B: INEQUALITIES [[0,1,0,0],[0,0,1,0]]
701: % $BF~NO$O(B polymake $B$N(B tree (doPolymake $B$N(B 1 get)
702: %>
703: /getLinearitySubspace {
704: /arg1 set
705: [/pdata /vv /ineq /rr /ii] pushVariables
706: [
707: /pdata arg1 def
708: {
709: /rr [ ] def
710: % POINTED $B$J$i(B max lin subspace $B$O(B 0.
711: pdata (POINTED) getNode tag 0 eq { } { exit} ifelse
712:
713: pdata (INEQUALITIES) getNode 2 get 0 get /ineq set
714: pdata (VERTICES) getNode 2 get 0 get /vv set
715: 0 1 vv length 1 sub {
716: /ii set
717: % -vv[ii] $B$,(B ineq $B$rK~$?$9$+D4$Y$k(B.
718: vv ii get ineq isInLinearSpace {
719: rr [vv ii get] join /rr set
720: } { } ifelse
721: } for
722: exit
723: } loop
724: /arg1 rr def
725: ] pop
726: popVariables
727: arg1
728: } def
729:
730: %<
731: % Usages: mm asir_matrix_image
732: % $B@8@.85$h$j@~7A6u4V$N4pDl$rF@$k(B.
733: %>
734: /asir_matrix_image {
735: /arg1 set
736: [/mm /rr] pushVariables
737: [(CurrentRingp)] pushEnv
738: [
739: /mm arg1 def
740: mm to_univNum /mm set
741: oxasir.ccc [ ] eq {
742: (Starting ox_asir server.) message
743: ox_asirConnectMethod
744: } { } ifelse
745: {
746: oxasir.ccc [(matrix_image) mm] asir
747: /rr set
748: rr null_to_zero /rr set
749: exit
750:
751: (asir_matrix_image: not implemented) error exit
752: } loop
753:
754: rr numerator /rr set
755: /arg1 rr def
756: ] pop
757: popEnv
758: popVariables
759: arg1
760: } def
761: [(asir_matrix_image)
762: [(Calling the function matrix_image of asir. It gets a reduced basis of a given matrix.)
763: (Example: [[1 2 3] [2 4 6]] asir_matrix_image)
764: ]] putUsages
765:
766: %<
767: % Usages: mm asir_matrix_kernel
768: % $BD>8r$9$k6u4V$N4pDl(B.
769: %>
770: /asir_matrix_kernel {
771: /arg1 set
772: [/mm /rr] pushVariables
773: [(CurrentRingp)] pushEnv
774: [
775: /mm arg1 def
776: mm to_univNum /mm set
777: oxasir.ccc [ ] eq {
778: (Starting ox_asir server.) message
779: ox_asirConnectMethod
780: } { } ifelse
781: {
782: oxasir.ccc [(matrix_kernel) mm] asir
783: /rr set
784: rr null_to_zero /rr set
785: exit
786:
787: (asir_matrix_image: not implemented) error exit
788: } loop
789: rr 1 get numerator /rr set
790: /arg1 rr def
791: ] pop
792: popEnv
793: popVariables
794: arg1
795: } def
796: [(asir_matrix_kernel)
797: [(Calling the function matrix_kernel of asir.)
798: (It gets a reduced basis of the kernel of a given matrix.)
799: (Example: [[1 2 3] [2 4 6]] asir_matrix_kernel)
800: ]] putUsages
801:
802: %<
803: % Usages: v null_to_zero
804: %>
805: /null_to_zero {
806: /arg1 set
807: [/pp /rr] pushVariables
808: [
809: /pp arg1 def
810: {
811: /rr pp def
812: pp isArray {
813: pp {null_to_zero} map /rr set
814: exit
815: }{ } ifelse
816:
817: pp tag 0 eq {
818: /rr (0).. def
819: exit
820: }{ } ifelse
821: exit
822: } loop
823: /arg1 rr def
824: ] pop
825: popVariables
826: arg1
827: } def
828: [(null_to_zero)
829: [(obj null_to_zero rob)
830: $It translates null to (0)..$
831: ]] putUsages
832:
1.4 takayama 833: %<
834: % Usages: newVector.with-1
835: % (-1).. $B$GKd$a$?%Y%/%H%k$r:n$k(B.
836: %>
837: /newVector.with-1 {
838: newVector { pop (-1).. } map
839: } def
840:
841:
1.1 takayama 842: % [2 0] lcm $B$O(B 0 $B$r$b$I$9$,$$$$$+(B? --> OK.
843:
844: %<
845: % Usages: mm addZeroForPolymake
846: % $B0J2<$NFs$D$N4X?t$O(B, toQuotientSpace $B$K$bMxMQ(B.
847: % Polymake INEQUALITIES $BMQ$K(B 0 $B$r;O$a$KB-$9(B.
848: % $BF~NO$O(B $B%j%9%H$N%j%9%H(B
849: % [[1,2], [3,4],[5,6]] --> [[0,1,2],[0,3,4],[0,5,6]]
850: %>
851: /addZeroForPolymake {
852: /arg1 set
853: [/mm /rr] pushVariables
854: [
855: /mm arg1 def
856: mm to_univNum /mm set
857: mm { [(0)..] 2 1 roll join } map /mm set
858: /arg1 mm def
859: ] pop
860: popVariables
861: arg1
862: } def
863:
864: %<
865: % Usages: mm cone.appendZero
866: %>
867: /cone.appendZero {
868: /arg1 set
869: [/mm /rr] pushVariables
870: [
871: /mm arg1 def
872: mm to_univNum /mm set
873: mm { [(0)..] join } map /mm set
874: /arg1 mm def
875: ] pop
876: popVariables
877: arg1
878: } def
879:
880: %<
881: % Usages: mm removeFirstFromPolymake
882: % $B;O$a$N(B 0 $B$r<h$j=|$/(B.
883: % $BF~NO$O(B $B%j%9%H$N%j%9%H(B
884: % [[0,1,2],[0,3,4],[0,5,6]] ---> [[1,2], [3,4],[5,6]]
885: %>
886: /removeFirstFromPolymake {
887: /arg1 set
888: [/mm /rr] pushVariables
889: [
890: /mm arg1 def
891: mm to_univNum /mm set
892: mm {rest} map /mm set
893: /arg1 mm def
894: ] pop
895: popVariables
896: arg1
897: } def
898:
899: %<
900: % Usages: mm genUnit
901: % [1,0,0,...] $B$r2C$($k$?$a$K@8@.(B.
902: % [[0,1,2], [0,3,4],[0,5,6]]--> [1,0,0]
903: %>
904: /genUnit {
905: /arg1 set
906: [/mm /rr /i] pushVariables
907: [
908: /mm arg1 def
909: mm 0 get length newVector /rr set
910: rr null_to_zero /rr set
911: rr 0 (1).. put
912: /arg1 rr def
913: ] pop
914: popVariables
915: arg1
916: } def
917:
918: %<
919: % Usages: mm genUnitMatrix
920: % [[0,1,2], [0,3,4],[0,5,6]]--> [[1,0,0],[0,1,0],[0,0,1]]
921: %>
922: /genUnitMatrix {
923: /arg1 set
924: [/mm /rr /nn /i] pushVariables
925: [
926: /mm arg1 def
927: mm 0 get length /nn set
928: [
929: 0 1 nn 1 sub {
930: /i set
931: nn newVector null_to_zero /mm set
932: mm i (1).. put
933: mm
934: } for
935: ]
936: /arg1 set
937: ] pop
938: popVariables
939: arg1
940: } def
941:
942: %<
943: %%note: 2004, 8/29 (sun)
944: % toQuotientSpace : Linearity space $B$G3d$k(B.
945: % Usages: ineq mm toQuotientSpace
946: % $BF~NO$O(B coneEq $B$N=PNO(B ineq
947: % $B$*$h$S(B doPolymake --> getLinearitySubspace ==> L
948: % [L,[1,0,0,...]] asir_matrix_kernel removeFirstFromPolymake $B$GF@$i$l$?(B mm
949: % $B=PNO$+$i(B 0 $B%Y%/%H%k$O:o=|(B.
950: % $B=PNO$b(B coneEq $B7A<0(B. $BFC$K(B polymake $BMQ$K(B 0 $B$r2C$($k$N$,I,MW(B.
951: % ref: getUnit, removeFirstFromPolymake, addZeroForPolymake,
952: % asir_matrix_kernel, getLinearitySubspace
953: %>
954: /toQuotientSpace {
955: /arg2 set
956: /arg1 set
957: [/ineq /mm /rr] pushVariables
958: [
959: /ineq arg1 def
960: /mm arg2 def
961:
962: ineq mm transpose mul /rr set
963:
964: /arg1 rr def
965: ] pop
966: popVariables
967: arg1
968: } def
969:
970: /test5.data
971: $polymake.data(polymake.INEQUALITIES([[0,1,-1,1,-1,0],[0,0,-1,0,-1,2],[0,0,-1,0,-1,2],[0,0,-2,0,-2,4],[0,-1,0,-1,0,2],[0,-2,0,-2,0,4]]),polymake.VERTICES([[0,0,-1,0,0,0],[0,-1,-1,0,0,0],[0,1,0,-1,0,0],[0,-1,0,1,0,0],[0,0,1,0,-1,0],[0,0,-1,0,1,0],[0,-2,-2,0,0,-1],[0,2,2,0,0,1]]),polymake.FACETS([[0,1,-1,1,-1,0],[0,-1,0,-1,0,2]]),polymake.AFFINE_HULL(),polymake.FEASIBLE(),polymake.NOT__POINTED(),polymake.FAR_FACE([polymake._set([0,1,2,3,4,5,6,7])]),polymake.VERTICES_IN_INEQUALITIES([polymake._set([1,2,3,4,5,6,7]),polymake._set([2,3,4,5,6,7]),polymake._set([2,3,4,5,6,7]),polymake._set([2,3,4,5,6,7]),polymake._set([0,2,3,4,5,6,7]),polymake._set([0,2,3,4,5,6,7])]),polymake.DIM([[5]]),polymake.AMBIENT_DIM([[5]]))$
972: def
973: %<
974: % Usages: test5
975: %% getConeInfo $B$rJQ99$9$l$P(B polymake $B$r8F$P$:$K%F%9%H$G$-$k(B.
976: %>
977: /test5 {
978: % test3b $B$h$j(B
979: /ww [(Dx) 1 (Dy) 2] def
980: % /ww [(x) 1 (y) -2 (Dx) 3 (Dy) 6] def
981: [(x,y) ring_of_differential_operators
982: [ww] weight_vector
983: 0] define_ring
984: [ (x Dx + y Dy -1).
985: (y^2 Dy^2 + 2 + y Dy ).
986: ] /gg set
987: gg {homogenize} map /gg set
988: [(AutoReduce) 1] system_variable
989: [gg] groebner 0 get /gg set
990: ww message
991:
992: ww gg coneEq getConeInfo /rr set
993: (Type in rr 0 get :: ) message
994: } def
995: %[5, [[1,0,1,0,-2],[0,1,0,1,-2]], $NOT__POINTED$ ]
996: % $B$3$N>l9g$O(B 2 $B<!85$^$GMn$9$H(B pointed cone $B$K$J$k(B.
997: % coneEq mmc transpose $B$r$b$H$K(B FACETS $B$r7W;;$9$l$P$h$$(B.
998:
999: %<
1000: % Usage: ceq getConeInfo
1001: % vw $B$O(B [v1 w1 v2 w2 ... ] $B7A<0(B. (v-w $B7A<0(B) w1, w2 $B$O(B univNumber $B$G$b$$$$(B.
1002: % g $B$O(B reduced Grobner basis $B$H$7$F(B vw g coneEq $B$r7W;;(B. $B$3$l$r(B getConeInfo $B$X(B.
1003: % Grobner cone $B$N(B $B<!85(B cdim (DIM), $BJd6u4V(B (linearity space ) $B$X$N9TNs(B mmc
1004: % linearity space $B<+BN(B, pointed or not__pointed
1005: % $B$D$^$j(B [cdim, L', L, PointedQ]
1006: % $B$r7W;;$7$FLa$9(B. (polymake $B7A<0$NM>J,$JItJ,$J$7(B)
1007: % polymake $BI,MW(B.
1008: % ref: coneEq
1009: % Global:
1010: % cone.getConeInfo.rr0, cone.getConeInfo.rr1 $B$K(B polymake $B$h$j$NLa$jCM$,$O$$$k(B.
1011: %>
1012: /getConeInfo {
1013: /arg1 set
1014: [/ww /g /ceq /ceq2 /cdim /mmc /mmL /rr /ineq /ppt] pushVariables
1015: [
1016: /ceq arg1 def
1017: ceq pruneZeroVector /ceq set
1018: ceq genPo2 /ceq2 set
1019: % ceq2 $B$O(B polymake.data(polymake.INEQUALITIES(...)) $B7A<0(B
1020: % polymake $B$G(B ceq2 $B$N<!85$N7W;;(B.
1021: /getConeInfo.ceq ceq def /getConeInfo.ceq2 ceq2 def
1022:
1023: cone.debug { (Calling polymake DIM.) message } { } ifelse
1024: [(DIM) ceq2] doPolymake 1 get /rr set
1025: cone.debug {(Done.) message } { } ifelse
1026: % test5 $B$K$O<!$N%3%a%s%H$H$j$5$k(B. $B>e$N9T$r%3%a%s%H%"%&%H(B.
1027: % test5.data tfbToTree /rr set
1028: /cone.getConeInfo.rr0 rr def
1029:
1030: rr (DIM) getNode /cdim set
1031: cdim 2 get 0 get 0 get 0 get to_univNum /cdim set
1032: % polymake $B$N(B DIM $B$O0l$D>.$5$$$N$G(B 1 $BB-$9(B.
1033: cdim (1).. add /cdim set
1034:
1035: rr (FACETS) getNode tag 0 eq {
1036: % FACETS $B$r;}$C$F$$$J$$$J$i:FEY7W;;$9$k(B.
1037: % POINTED, NOT__POINTED $B$bF@$i$l$k(B
1038: cone.debug { (Calling polymake FACETS.) message } { } ifelse
1039: [(FACETS) ceq2] doPolymake 1 get /rr set
1040: cone.debug { (Done.) message } { } ifelse
1041: } { } ifelse
1042:
1043: rr (VERTICES) getNode tag 0 eq {
1044: (internal error: VERTICES is not found.) error
1045: } { } ifelse
1046:
1047: /cone.getConeInfo.rr1 rr def
1048:
1049: rr (NOT__POINTED) getNode tag 0 eq {
1050: % cone $B$,(B pointed $B$N;~$O(B mmc $B$OC10L9TNs(B. genUnitMatrix $B$r;H$&(B.
1051: % VERTICES $B$h$j0l$D>.$5$$%5%$%:(B.
1052: /mmc
1053: [ rr (VERTICES) getNode 2 get 0 get 0 get rest]
1054: genUnitMatrix
1055: def
1056: /mmL [ ] def
1057: /ppt (POINTED) def
1058: } {
1059: % pointed $B$G$J$$>l9g(B,
1060: % cone $B$N@~7AItJ,6u4V$r7W;;(B.
1061: rr getLinearitySubspace /mmL set
1062: [mmL genUnit] mmL join /mmc set % [1,0,0,...] $B$rB-$9(B.
1063: mmc asir_matrix_kernel /mmc set % $BJd6u4V(B
1064: mmc removeFirstFromPolymake /mmc set % $B$R$H$D>.$5$$%5%$%:$K(B.
1065:
1066: [mmL genUnit] mmL join asir_matrix_image
1067: removeFirstFromPolymake /mmL set
1068: mmL asir_matrix_image /mmL set % Linearity space $B$r5a$a$k(B. rm 0vector
1069: /ppt (NOT__POINTED) def
1070: } ifelse
1071: /arg1 [[cdim mmc mmL ppt] rr] def
1072: ] pop
1073: popVariables
1074: arg1
1075: } def
1076:
1077:
1078: /test.put {
1079: /dog [(dog) [[(legs) 4] ] [1 2 3 ]] [(class) (tree)] dc def
1080: /man [(man) [[(legs) 2] ] [1 2 3 ]] [(class) (tree)] dc def
1081: /ma [(mammal) [ ] [man dog]] [(class) (tree)] dc def
1082: /fan [ma 1 copy] def
1083: ma (dog) getNode /dd set
1084: dd 2 get /dd2 set
1085: dd2 1 0 put
1086: ma message
1087:
1088: fan message
1089: } def
1090:
1091: /test6.data
1092: $polymake.data(polymake.INEQUALITIES([[0,1,-1,1,-1,0],[0,0,-1,0,-1,2],[0,0,-1,0,-1,2],[0,0,-2,0,-2,4],[0,-1,0,-1,0,2],[0,-2,0,-2,0,4]]),polymake.VERTICES([[0,0,-1,0,0,0],[0,-1,-1,0,0,0],[0,1,0,-1,0,0],[0,-1,0,1,0,0],[0,0,1,0,-1,0],[0,0,-1,0,1,0],[0,-2,-2,0,0,-1],[0,2,2,0,0,1]]),polymake.FACETS([[0,1,-1,1,-1,0],[0,-1,0,-1,0,2]]),polymake.AFFINE_HULL(),polymake.FEASIBLE(),polymake.NOT__POINTED(),polymake.FAR_FACE([polymake._set([0,1,2,3,4,5,6,7])]),polymake.VERTICES_IN_INEQUALITIES([polymake._set([1,2,3,4,5,6,7]),polymake._set([2,3,4,5,6,7]),polymake._set([2,3,4,5,6,7]),polymake._set([2,3,4,5,6,7]),polymake._set([0,2,3,4,5,6,7]),polymake._set([0,2,3,4,5,6,7])]))$
1093: def
1094: % tfbToTree
1095:
1096: /arrayToTree { [(class) (tree)] dc } def
1097:
1098: %<
1099: % polymake $B$h$jF@$i$l$?(B TreeObject $B$+$i(B TreeObject cone $B$r@8@.$9$k(B.
1100: % Usages: test6.data tfbToTree newCone $B$GF0:n%F%9%H(B
1101: %>
1102: /test6 {
1103: test6.data tfbToTree /rr set
1104: rr newCone /rr2 set
1105: } def
1106:
1107: %<
1108: % Usages: doPolymakeObj newCone
1109: %>
1110: /newCone {
1111: /arg1 set
1112: [/polydata /cone /facets /vertices /flipped /ineq
1113: /facetsv /rr] pushVariables
1114: [
1115: /polydata arg1 def
1116: polydata (FACETS) getNode tag 0 eq {
1117: (newCone : no FACETS data.) error
1118: } { } ifelse
1119: % facets $B$OM-M}?t$N>l9g@55,2=$9$k(B. data/test11 $B$G(B $BM-M}?t$G$k(B.
1120: polydata (FACETS) getNode 2 get 0 get to_univNum
1121: { nnormalize_vec} map /facets set
1122: [[ ] ] facets join shell rest removeFirstFromPolymake /facets set
1.2 takayama 1123: facets length 0 eq
1124: {(Internal error. Facet data is not obtained. See OpenXM_tmp.) error} { } ifelse
1.1 takayama 1125: % vertices $B$O(B cone $B$N>e$K$"$k$N$G@0?tG\(B OK. $B@55,$+$9$k(B.
1126: polydata (VERTICES) getNode 2 get 0 get to_univNum
1127: { nnormalize_vec} map /vertices set
1128: [[ ] ] vertices join shell rest removeFirstFromPolymake /vertices set
1129: % inequalities $B$OM-M}?t$N>l9g@55,2=$9$k(B.
1130: polydata (INEQUALITIES) getNode 2 get 0 get to_univNum
1131: { nnormalize_vec } map /ineq set
1132: [[ ] ] ineq join shell rest removeFirstFromPolymake /ineq set
1133:
1.4 takayama 1134: % nextcid, nextfid $B$r2C$($k(B. nextcid $B$O(B nextConeId $B$NN,(B. $B$H$J$j$N(B cone $BHV9f(B.
1135: % nextfid $B$O(B nextFacetId $B$NN,(B. $B$H$J$j$N(B cone $B$N(B facet
1136: % $BHV9f(B.
1.1 takayama 1137: [(cone) [ ]
1138: [
1139: [(facets) [ ] facets] arrayToTree
1140: [(flipped) [ ] facets length newVector null_to_zero] arrayToTree
1141: [(facetsv) [ ] facets vertices newCone_facetsv] arrayToTree
1.4 takayama 1142: [(nextcid) [ ] facets length newVector.with-1 ] arrayToTree
1143: [(nextfid) [ ] facets length newVector.with-1 ] arrayToTree
1.1 takayama 1144: [(vertices) [ ] vertices] arrayToTree
1145: [(inequalities) [ ] ineq] arrayToTree
1146: ]
1147: ] arrayToTree /cone set
1148: /arg1 cone def
1149: ] pop
1150: popVariables
1151: arg1
1152: } def
1153:
1154: %<
1155: % Usages: newCone_facetv
1156: % facet vertices newCone_facetv
1157: % facet $B$K$N$C$F$$$k(B vertices $B$r$9$Y$FNs5s(B.
1158: %>
1159: /newCone_facetv {
1160: /arg2 set
1161: /arg1 set
1162: [/facet /vertices] pushVariables
1163: [
1164: /facet arg1 def /vertices arg2 def
1165: [
1166: 0 1 vertices length 1 sub {
1167: /ii set
1168: facet vertices ii get mul isZero
1169: { vertices ii get } { } ifelse
1170: } for
1171: ]
1172: /arg1 set
1173: ] pop
1174: popVariables
1175: arg1
1176: } def
1177:
1178: %<
1179: % Usages: newCone_facetsv
1180: % facets vertices newCone_facetv
1181: % facets $B$K$N$C$F$$$k(B vertices $B$r$9$Y$FNs5s(B. $B%j%9%H$r:n$k(B.
1182: %>
1183: /newCone_facetsv {
1184: /arg2 set
1185: /arg1 set
1186: [/facets /vertices] pushVariables
1187: [
1188: /facets arg1 def /vertices arg2 def
1189: facets { vertices newCone_facetv } map
1190: /arg1 set
1191: ] pop
1192: popVariables
1193: arg1
1194: } def
1195:
1196: %<
1.2 takayama 1197: % Usages: [gb weight] newConeGB
1198: % gb $B$H(B weight $B$r(B tree $B7A<0$K$7$F3JG<$9$k(B.
1199: %>
1200: /newConeGB {
1201: /arg1 set
1202: [/gbdata /gg /ww /rr] pushVariables
1203: [
1204: /gbdata arg1 def
1205: % gb
1206: gbdata 0 get /gg set
1207: % weight
1208: gbdata 1 get /ww set
1209: %
1210: [(coneGB) [ ]
1211: [
1212: [(grobnerBasis) [ ] gg] arrayToTree
1213: [(weight) [ ] [ww]] arrayToTree
1214: [(initial) [ ] gg { ww 2 get weightv init } map ] arrayToTree
1215: ]
1216: ] arrayToTree /rr set
1217: /arg1 rr def
1218: ] pop
1219: popVariables
1220: arg1
1221: } def
1222:
1223: %<
1.1 takayama 1224: % Usages: cone_random
1225: %>
1226: /cone_random.start (2).. def
1227: /cone_random {
1228: [(tdiv_qr)
1229: cone_random.start (1103515245).. mul
1230: (12345).. add
1231:
1232: (2147483646)..
1233: ] mpzext 1 get /cone_random.start set
1234: cone_random.start
1235: } def
1236:
1237: /cone_random.limit 40 def
1238: /cone_random_vec {
1239: /arg1 set
1240: [/nn /rr] pushVariables
1241: [
1242: /nn arg1 def
1243: [
1244: 0 1 nn 1 sub {
1245: pop
1246: [(tdiv_qr) cone_random cone_random.limit] mpzext 1 get
1247: } for
1248: ] /arg1 set
1249: ] pop
1250: popVariables
1251: arg1
1252: } def
1253:
1254: %<
1255: % Usages: getNewRandomWeight
1256: %% max dim $B$N(B cone $B$r@8@.$9$k$?$a$K(B, random $B$J(B weight $B$r@8@.$9$k(B.
1257: %% h, H $B$N=hM}$bI,MW(B.
1258: %% $B@)Ls>r7o(B u+v >= 2t $B$r$_$?$9(B weight $B$,I,MW(B. $B$3$l$r$I$N$h$&$K:n$k$N$+(B?
1259: %>
1260: /getNewRandomWeight {
1261: /arg1 set
1262: [/vv /vvd /rr] pushVariables
1263: [
1264: /vv arg1 def
1265: vv { (D) 2 1 roll 2 cat_n } map /vvd set
1266: ] pop
1267: popVariables
1268: arg1
1269: } def
1270:
1271: % test7 : univNum $B$N(B weight $B$,@5$7$/G'<1$5$l$k$+$N%F%9%H(B
1272: % aux-cone.sm1
1273:
1274: %<
1275: % Usages: n d coneEqForSmallFan.2 (cone.type 2 $B@lMQ(B: x,y,Dx,Dy,h)
1276: % n $BJQ?t$N?t(B, d zero $B$K$7$J$$JQ?t$N?t(B. d $B$O(B max dim cone $B$N<!85$H$J$k(B.
1277: % $B$O$8$a$+$i(B d $B8D$NJQ?t(B.
1278: % 4, 2 , s,t,x,y $B$J$i(B weight $B$O(B s,t,Ds,Dt $B$N$_(B.
1279: % u_i + v_i >= 0 , u_i = v_i = 0.
1280: % homog $BJQ?t$N>r7o(B u_i+v_i >= t, i.e, -t >= 0 $B$bF~$l$k(B.
1281: % coneEq $B$N7k2L$H(B coneEqForSmallFan.2 $B$N7k2L$r(B join $B$7$F(B
1282: % getConeInfo or newCone
1283: % note-cone.sm1 2004.8.31 $B$r8+$h(B. w_ineq $B$"$?$j(B.
1284: % cone.local $B$,@_Dj$5$l$F$$$k$H(B u_i <= 0 $B$b>r7o$KF~$k(B.
1285: %>
1286: /coneEqForSmallFan.2 {
1287: /arg2 set
1288: /arg1 set
1289: [/n /d /nn /dd /ii /tt] pushVariables
1290: [
1291: /n arg1 def
1292: /d arg2 def
1293: n to_int32 /n set
1294: d to_int32 /d set
1295: /dd n d add def
1296: /nn n n add def
1297:
1298: % 0 ~ d-1, n ~ dd-1 $B$G$O(B u_i + v_i = 0
1299: % d ~ n-1, dd ~ nn-1 $B$G$O(B u_i=v+i = 0.
1300: % -t >= 0
1301: [
1302: % d ~ n-1, dd ~ nn-1 $B$G$O(B u_i=v+i = 0.
1303: d 1 n 1 sub {
1304: /ii set
1305: % [ 0,0, ..., 0,1,0,... ; 0] $B$r@8@.(B
1306: nn 1 add newVector null_to_zero /tt set
1307: tt ii (1).. put
1308: tt
1309: % [ 0,0, ..., 0,-1,0,... ; 0] $B$r@8@.(B
1310: nn 1 add newVector null_to_zero /tt set
1311: tt ii (-1).. put
1312: tt
1313: } for
1314: dd 1 nn 1 sub {
1315: /ii set
1316: nn 1 add newVector null_to_zero /tt set
1317: tt ii (1).. put
1318: tt
1319: nn 1 add newVector null_to_zero /tt set
1320: tt ii (-1).. put
1321: tt
1322: } for
1323:
1324: % 0 ~ d-1, n ~ dd-1 $B$G$O(B u_i + v_i = 0
1325: 0 1 d 1 sub {
1326: /ii set
1327: nn 1 add newVector null_to_zero /tt set
1328: tt ii (1).. put
1329: tt ii n add (1).. put
1330: tt
1331:
1332: nn 1 add newVector null_to_zero /tt set
1333: tt ii (-1).. put
1334: tt ii n add (-1).. put
1335: tt
1336:
1337: } for
1338:
1339: % -t >= 0
1340: cone.h0 {
1341: % t = 0
1342: nn 1 add newVector null_to_zero /tt set
1343: tt nn (1).. put
1344: tt
1345: nn 1 add newVector null_to_zero /tt set
1346: tt nn (-1).. put
1347: tt
1348: }
1349: {
1350: % -t >= 0
1351: nn 1 add newVector null_to_zero /tt set
1352: tt nn (-1).. put
1353: tt
1354: } ifelse
1355:
1356: % cone.local $B$,(B 1 $B$N;~(B
1357: % 0 ~ d-1 $B$G$O(B -u_i >= 0
1358: cone.local {
1359: 0 1 d 1 sub {
1360: /ii set
1361: nn 1 add newVector null_to_zero /tt set
1362: tt ii (-1).. put
1363: tt
1364: } for
1365: } { } ifelse
1366: ] /rr set
1367: /arg1 rr to_univNum def
1368: ] pop
1369: popVariables
1370: arg1
1371: } def
1372:
1373: %<
1374: % Usages: n d coneEqForSmallFan.1 (cone.type 1 $B@lMQ(B: x,y,Dx,Dy,h,H)
1375: % cone.type 2 $B$G$O(B x,y,Dx,Dy,h
1376: % coneEqForSmallFan.2 $B$N7k2L$rMQ$$$F@8@.(B.
1377: % H $B$N>r7o$r2C$($k(B.
1378: %>
1379: /coneEqForSmallFan.1 {
1380: /arg2 set
1381: /arg1 set
1382: [/n /d /i /j /rr /tt /tt2] pushVariables
1383: [
1384: /n arg1 def /d arg2 def
1385: n d coneEqForSmallFan.2 /rr set
1386: rr cone.appendZero /rr set
1387: % H $BMQ$N(B 0 $B$r2C$($k(B.
1388: % $B$H$j$"$($:(B t' = 0 $B$G$-$a$&$A(B.
1389: cone.h0 { } { (cone.h0 = 0 has not yet been implemented.) error } ifelse
1390: n 2 mul 2 add newVector null_to_zero /tt set
1391: tt n 2 mul 2 add 1 sub (-1).. put
1392: n 2 mul 2 add newVector null_to_zero /tt2 set
1393: tt2 n 2 mul 2 add 1 sub (1).. put
1394: rr [tt tt2] join /rr set
1395: /arg1 rr to_univNum def
1396: ] pop
1397: popVariables
1398: arg1
1399: } def
1400:
1401: %<
1402: % Usages: vv ineq toQuotientCone
1403: % weight space $B$N(B $B%Q%i%a!<%?$D$1$N$?$a$K;H$&(B.
1404: % cone.V $B$r5a$a$?$$(B. vv $B$O(B doPolymakeObj (VERTICES) getNode 2 get 0 get $B$GF@$k(B.
1405: % vertices $B$N(B non-negative combination $B$,(B cone.
1406: % vertice cone.w_ineq isInLinearSubspace $B$J$i<h$j=|$/(B.
1407: % $B$D$^$j(B vertice*cone.w_ineq = 0 $B$J$i<h$j=|$/(B.
1408: %
1409: % $B$3$l$G@5$7$$(B? $B>ZL@$O(B? $B$^$@ESCf(B. cone.W $B$r5a$a$k$N$K;H$&(B. (BUG)
1410: % cone.w_cone 1 get (VERTICES) getNode :: $B$HHf3S$;$h(B.
1411: % $B$3$N4X?t$r8F$s$G(B cone.W $B$r:n$k$N$OITMW$+$b(B.
1412: %
1413: % Example: cf. parametrizeSmallFan
1414: % 4 2 coneEqForSmallFan.2 /cone.w_ineq set cone.w_ineq getConeInfo /rr set
1415: % rr 1 get (VERTICES) getNode 2 get 0 get removeFirstFromPolymake /vv set
1416: % vv cone.w_ineq toQuotientCone pmat
1417: %>
1418: /toQuotientCone {
1419: /arg2 set /arg1 set
1420: [/vv /ineq /rr] pushVariables
1421: [
1422: /vv arg1 def /ineq arg2 def
1423: vv {
1424: dup
1425: ineq isInLinearSpace 1 eq { pop }
1426: { } ifelse
1427: } map /arg1 set
1428: ] pop
1429: popVariables
1430: arg1
1431: } def
1432:
1433: %<
1434: % Usages: n d parametrizeSmallFan
1435: % n : x $BJQ?t$N?t(B.
1436: % d : 0 $B$K$7$J$$(B weight $B$N?t(B.
1437: % $B<!$NBg0hJQ?t$b@_Dj$5$l$k(B.
1438: % cone.W : weight $B$r%Q%i%a!<%?$E$1$9$k%Y%/%H%k$NAH(B.
1439: % cone.Wpos : i $B$,(B 0 ~ Wpos-1 $B$NHO0O$N$H$-(B V[i] $B$X$O(B N $B$N85$r3]$1;;$7$F$h$$(B,
1440: % i $B$,(B Wpos ~ $B$NHO0O$N$H$-(B V[i] $B$X$O(B Z $B$N85$r3]$1;;$7$F$h$$(B.
1441: % cone.w_ineq : weight space $B$NITEy<0@)Ls(B. $B0J8e$N7W;;$G>o$KIU2C$9$k(B.
1442: % cone.w_cone : w_ineq $B$r(B polymake $B$G(B getConeInfo $B$7$?7k2L(B.
1443: % Example: /cone.local 1 def ; 4 2 parametrizeSmallFan pmat
1444: % Example: /cone.local 0 def ; 4 2 parametrizeSmallFan pmat
1445: %>
1446: /parametrizeSmallFan {
1447: /arg2 set /arg1 set
1448: [/n /d /vv /coneray] pushVariables
1449: [
1450: /n arg1 def /d arg2 def
1451: {
1452: cone.type 1 eq {
1453: n d coneEqForSmallFan.1 /cone.w_ineq set
1454: exit
1455: } { } ifelse
1456: cone.type 2 eq {
1457: n d coneEqForSmallFan.2 /cone.w_ineq set
1458: exit
1459: } { } ifelse
1460: (This cone.type has not yet been implemented.) error
1461: } loop
1462: cone.w_ineq getConeInfo /cone.w_cone set
1463: cone.w_cone 1 get (VERTICES) getNode 2 get 0 get
1464: removeFirstFromPolymake /vv set
1465:
1466: vv cone.w_ineq toQuotientCone /coneray set
1467: coneray length /cone.Wpos set
1468:
1469: coneray cone.w_cone 0 get 2 get join /cone.W set
1470: /arg1 cone.W def
1471: ] pop
1472: popVariables
1473: arg1
1474: } def
1475:
1476: %<
1477: % Usages: n d coneEqForTotalFan.2 (cone.type 2 $B@lMQ(B: x,y,Dx,Dy,h)
1478: % n $BJQ?t$N?t(B,
1479: % d 0 $B$K$7$J$$JQ?t(B.
1480: % u_i + v_i >= 0 ,
1481: % homog $BJQ?t$N>r7o(B u_i+v_i >= 0, t = 0 $B$bF~$l$k(B.
1482: % coneEq $B$N7k2L$H(B coneEqForSmallFan.2 $B$N7k2L$r(B join $B$7$F(B
1483: % getConeInfo or newCone
1484: % cone.local $B$,@_Dj$5$l$F$$$k$H(B u_i <= 0 $B$b>r7o$KF~$k(B.
1485: %>
1486: /coneEqForTotalFan.2 {
1487: /arg2 set
1488: /arg1 set
1489: [/n /nn /dd /ii /tt] pushVariables
1490: [
1491: /n arg1 def
1492: /d arg2 def
1493: n to_int32 /n set
1494: d to_int32 /d set
1495: /nn n n add def
1496: /dd n d add def
1497:
1498: % 0 ~ d-1, n ~ dd-1 $B$G$O(B u_i + v_i >= 0
1499: % d ~ n-1, dd ~ nn-1 $B$G$O(B u_i=v+i = 0.
1500: % t = 0
1501: [
1502: % d ~ n-1, dd ~ nn-1 $B$G$O(B u_i=v+i = 0.
1503: d 1 n 1 sub {
1504: /ii set
1505: % [ 0,0, ..., 0,1,0,... ; 0] $B$r@8@.(B
1506: nn 1 add newVector null_to_zero /tt set
1507: tt ii (1).. put
1508: tt
1509: % [ 0,0, ..., 0,-1,0,... ; 0] $B$r@8@.(B
1510: nn 1 add newVector null_to_zero /tt set
1511: tt ii (-1).. put
1512: tt
1513: } for
1514: dd 1 nn 1 sub {
1515: /ii set
1516: nn 1 add newVector null_to_zero /tt set
1517: tt ii (1).. put
1518: tt
1519: nn 1 add newVector null_to_zero /tt set
1520: tt ii (-1).. put
1521: tt
1522: } for
1523:
1524: % 0 ~ d-1, n ~ dd-1 $B$G$O(B u_i + v_i >= 0
1525: 0 1 d 1 sub {
1526: /ii set
1527: nn 1 add newVector null_to_zero /tt set
1528: tt ii (1).. put
1529: tt ii n add (1).. put
1530: tt
1531:
1532: } for
1533:
1534: % t = 0
1535: cone.h0 {
1536: % t = 0
1537: nn 1 add newVector null_to_zero /tt set
1538: tt nn (1).. put
1539: tt
1540: nn 1 add newVector null_to_zero /tt set
1541: tt nn (-1).. put
1542: tt
1543: }
1544: {
1545: (coneForTotalFan.2. Not implemented.) error
1546: } ifelse
1547:
1548: % cone.local $B$,(B 1 $B$N;~(B
1549: % 0 ~ d-1 $B$G$O(B -u_i >= 0
1550: cone.local {
1551: 0 1 d 1 sub {
1552: /ii set
1553: nn 1 add newVector null_to_zero /tt set
1554: tt ii (-1).. put
1555: tt
1556: } for
1557: } { } ifelse
1558: ] /rr set
1559: /arg1 rr to_univNum def
1560: ] pop
1561: popVariables
1562: arg1
1563: } def
1564:
1565: %<
1566: % Usages: n d parametrizeTotalFan
1567: % n : x $BJQ?t$N?t(B.
1568: % d : 0 $B$K$7$J$$?t(B.
1569: % $B<!$NBg0hJQ?t$b@_Dj$5$l$k(B.
1570: % cone.W : weight $B$r%Q%i%a!<%?$E$1$9$k%Y%/%H%k$NAH(B.
1571: % cone.Wpos : i $B$,(B 0 ~ Wpos-1 $B$NHO0O$N$H$-(B V[i] $B$X$O(B N $B$N85$r3]$1;;$7$F$h$$(B,
1572: % i $B$,(B Wpos ~ $B$NHO0O$N$H$-(B V[i] $B$X$O(B Z $B$N85$r3]$1;;$7$F$h$$(B.
1573: % cone.w_ineq : weight space $B$NITEy<0@)Ls(B. $B0J8e$N7W;;$G>o$KIU2C$9$k(B.
1574: % cone.w_ineq $B$r(B getConeInfo $B$7$?7k2L$O(B cone.w_cone
1575: % Example: /cone.local 1 def ; 3 parametrizeSmallFan pmat
1576: % Example: /cone.local 0 def ; 3 parametrizeSmallFan pmat
1577: % local $B$,(B 1 $B$@$H(B u_i <= 0 $B$K$J$k(B.
1578: %>
1579: /parametrizeTotalFan {
1580: /arg2 set
1581: /arg1 set
1582: [/n /d /vv /coneray] pushVariables
1583: [
1584: /n arg1 def /d arg2 def
1585: {
1586: cone.type 2 eq { n d coneEqForTotalFan.2 /cone.w_ineq set exit}
1587: { } ifelse
1588: (This cone.type has not yet been implemented.) error
1589: } loop
1590: cone.w_ineq getConeInfo /cone.w_cone set
1591: cone.w_cone 1 get (VERTICES) getNode 2 get 0 get
1592: removeFirstFromPolymake /vv set
1593:
1594: vv cone.w_ineq toQuotientCone /coneray set
1595: coneray length /cone.Wpos set
1596:
1597: coneray cone.w_cone 0 get 2 get join /cone.W set
1598: /arg1 cone.W def
1599: ] pop
1600: popVariables
1601: arg1
1602: } def
1603:
1604: %<
1605: % Usages: vlist wlist cone_wtowv
1606: % [x y Dx Dy h] [-1 0 1 0 0] ==> [(x) -1 (Dx) 1] $B$r:n$k(B.
1607: %>
1608: /cone_wtowv {
1609: /arg2 set /arg1 set
1610: [/vlist /wlist /ii] pushVariables
1611: [
1612: /vlist arg1 def
1613: /wlist arg2 def
1614: wlist length vlist length eq {
1615: } { (cone_wtowv: length of the argument must be the same.) error} ifelse
1616:
1617: wlist to_int32 /wlist set
1618: [
1619: 0 1 wlist length 1 sub {
1620: /ii set
1621: wlist ii get 0 eq { }
1622: { vlist ii get wlist ii get } ifelse
1623: } for
1624: ] /arg1 set
1625: ] pop
1626: popVariables
1627: arg1
1628: } def
1629:
1630: %<
1631: % Usages: pruneZeroVector
1632: % genPo, getConeInfo $BEy$NA0$K;H$&(B. 0 $B%Y%/%H%k$O0UL#$N$J$$@)Ls$J$N$G=|$/(B.
1.2 takayama 1633: % $BF1$8@)Ls>r7o$b$N$>$/(B. polymake FACET $B$,@5$7$/F0$+$J$$>l9g$,$"$k$N$G(B.
1634: % cf. pear/OpenXM_tmp/x3y2.poly, x^3+y^2, x^2+y^3 data/test15.sm1
1.1 takayama 1635: %>
1636: /pruneZeroVector {
1637: /arg1 set
1638: [/mm /ii /jj /tt] pushVariables
1639: [
1640: /mm arg1 def
1641: mm to_univNum /mm set
1.2 takayama 1642: [ [ ] ] mm join shell rest uniq /mm set
1.1 takayama 1643: [
1644: 0 1 mm length 1 sub {
1645: /ii set
1646: mm ii get /tt set
1647: {
1648: 0 1 tt length 1 sub {
1649: /jj set
1650: tt jj get (0).. eq { }
1651: { tt exit } ifelse
1652: } for
1653: exit
1654: } loop
1655: } for
1656: ] /arg1 set
1657: ] pop
1658: arg1
1659: } def
1660:
1661: %<
1662: % Usages: a projectIneq v , dim(a) = n, dim(v) = d
1663: % a*cone.Wt*cone.Lpt
1664: %>
1665: /projectIneq {
1666: cone.Wt mul cone.Lpt mul
1667: } def
1668:
1669: %<
1670: % Usages: v liftWeight [w vw], dim(v) = d, dim(w) = n, vw : vw $B7A<0$N(B weight
1671: % v*cone.Lp*cone.W cone.vlist w cone_wtowv
1672: %>
1673: /liftWeight {
1674: /arg1 set
1675: [/v /w /vw] pushVariables
1676: [
1677: /v arg1 def
1678: v cone.Lp mul cone.W mul /w set
1679: [w cone.vlist w cone_wtowv] /arg1 set
1680: ] pop
1681: popVariables
1682: arg1
1683: } def
1684:
1685: %<
1686: % Usage: m isZero
1687: % dr.sm1 $B$X0\$9(B.
1688: %>
1689: /isZero {
1690: /arg1 set
1691: [/mm /ans /ii] pushVariables
1692: [
1693: /mm arg1 def
1694: /ans 1 def
1695: mm isArray {
1696: 0 1 mm length 1 sub {
1697: /ii set
1698: mm ii get isZero /ans set
1699: ans 0 eq { exit } { } ifelse
1700: } for
1701: } {
1702: {
1703: mm tag 1 eq {/ans mm 0 eq def exit} { } ifelse
1704: mm isPolynomial { /ans mm (0). eq def exit } { } ifelse
1705: mm isUniversalNumber { /ans mm (0).. eq def exit } { } ifelse
1706: /ans 0 def exit
1707: } loop
1708: } ifelse
1709: /arg1 ans def
1710: ] pop
1711: popVariables
1712: arg1
1713: } def
1714: [(isZero)
1715: [(m isZero bool)]] putUsages
1716:
1717: %<
1718: % Usage: m isNonNegative
1719: % dr.sm1 $B$X0\$9(B.
1720: %>
1721: /isNonNegative {
1722: /arg1 set
1723: [/mm /ans /ii] pushVariables
1724: [
1725: /mm arg1 def
1726: /ans 1 def
1727: mm isArray {
1728: 0 1 mm length 1 sub {
1729: /ii set
1730: mm ii get isNonNegative /ans set
1731: ans 0 eq { exit } { } ifelse
1732: } for
1733: } {
1734: {
1735: mm tag 1 eq {/ans mm 0 gt mm 0 eq or def exit} { } ifelse
1736: mm isUniversalNumber { /ans mm (0).. gt mm (0).. eq or def exit }
1737: { } ifelse
1738: mm isRational { mm (numerator) dc mm (denominator) dc mul /mm set
1739: /ans mm (0).. gt mm (0).. eq or def exit } { } ifelse
1740: /ans 0 def exit
1741: } loop
1742: } ifelse
1743: /arg1 ans def
1744: ] pop
1745: popVariables
1746: arg1
1747: } def
1748: [(isNonNegative)
1749: [(m isNonNegative bool)
1750: (In case of matrix, m[i,j] >= 0 must hold for all i,j.)
1751: ]] putUsages
1752:
1753: % Global variable: cone.weightBorder
1754: % /cone.weightBorder null def $BITMW$G$"$m$&(B. getStartingCone $B$G@_Dj$5$l$k(B.
1755:
1756: %<
1757: % Usages: cone i isOnWeigthBorder
1758: % cone $B$N(B i $BHVL\$N(B facet $B$,(B weight $B6u4V$N6-3&$K$"$k$+(B?
1759: % $BBg0hJQ?t(B cone.weightBorder $B$,@_Dj$5$l$F$k$3$H(B.
1760: % $B$3$NJQ?t$O(B cone $B$N(B facet $B%Y%/%H%k$N%j%9%H(B.
1761: % $B$3$NJQ?t$O(B setWeightBorder $B$G@_Dj(B
1762: % cone.weightBorder[0] or cone.weightBorder[1] or ...
1763: % /ccone cone.startingCone def ccone 0 isOnWeightBorder
1764: % ccone 1 isOnWeightBorder
1765: %>
1766: /isOnWeightBorder {
1767: /arg2 set /arg1 set
1768: [/cone /facet_i /i /j /vv /co /ans] pushVariables
1769: [
1770: /cone arg1 def /facet_i arg2 def
1771: facet_i to_int32 /facet_i set
1772: /ans 0 def
1773: cone (facetsv) getNode 2 get facet_i get /vv set % Facet $B$r(B vertex $BI=8=(B.
1774: {
1775: 0 1 cone.weightBorder length 1 sub {
1776: /i set
1777: cone.weightBorder i get /co set % co $B$K@)Ls>r7o(B
1778: vv cone.Lp mul % vv $B$r(B weight space $B$X(B lift.
1779: co mul isZero
1780: { /ans 1 def exit } { } ifelse
1781: } for
1782: exit
1783: } loop
1784: /arg1 ans def
1785: ] pop
1786: popVariables
1787: arg1
1788: } def
1789:
1790: %<
1791: % Usages: cone i markFlipped
1792: % cone $B$N(B i $BHVL\$N(B facet $B$K(B flipped $B$N0u$r$D$1$k(B. cone $B<+BN$,JQ99$5$l$k(B.
1793: % cone $B$O(B class-tree. Constructor $B$O(B newCone
1794: %>
1795: /markFlipped {
1796: /arg2 set /arg1 set
1797: [/cone /facet_i /vv] pushVariables
1798: [
1799: /cone arg1 def /facet_i arg2 def
1800: facet_i to_int32 /facet_i set
1801: cone (flipped) getNode 2 get /vv set
1802: vv facet_i (1).. put
1803: ] pop
1804: popVariables
1805: } def
1806:
1.4 takayama 1807: %<
1808: % Usages: cone i [cid fid] markNext
1809: % cone $B$N(B i $BHVL\$N(B facet $B$N$H$J$j$N(B cone id (cid) $B$H(B face id (fid) $B$r@_Dj$9$k(B.
1810: % cone $B$N(B nextcid[i] = cid; nextfid[i] = fid $B$H$J$k(B.
1811: % cone $B<+BN$,JQ99$5$l$k(B.
1812: % cone $B$O(B class-tree.
1813: %>
1814: /markNext {
1815: /arg3 set /arg2 set /arg1 set
1816: [/cone /facet_i /vv /nextid] pushVariables
1817: [
1818: /cone arg1 def /facet_i arg2 def /nextid arg3 def
1819: facet_i to_int32 /facet_i set
1820: cone (nextcid) getNode 2 get /vv set
1821: vv facet_i , nextid 0 get to_univNum , put
1822:
1823: cone (nextfid) getNode 2 get /vv set
1824: vv facet_i , nextid 1 get to_univNum , put
1825: ] pop
1826: popVariables
1827: } def
1828:
1.1 takayama 1829:
1830:
1831: %<
1832: % Usages: cone getNextFacet i
1833: % flipped $B$N(B mark $B$N$J$$(B facet $B$N(B index facet_i $B$rLa$9(B.
1834: % $B$=$l$,$J$$$H$-$O(B null
1835: %>
1836: /getNextFacet {
1837: /arg1 set
1838: [/cone /facet_i /vv /ii] pushVariables
1839: [
1840: /cone arg1 def
1841: /facet_i null def
1842: cone (flipped) getNode 2 get /vv set
1843: 0 1 vv length 1 sub {
1844: /ii set
1845: vv ii get to_int32 0 eq { /facet_i ii def exit }
1846: { } ifelse
1847: } for
1848: /arg1 facet_i def
1849: ] pop
1850: popVariables
1851: arg1
1852: } def
1853:
1854: %<
1855: % Usages: cone i epsilon flipWeight
1856: % cone $B$N(B i $BHVL\$N(B facet $B$K$+$s$7$F(B flip $B$9$k(B.
1857: % $B?7$7$$(B weight $B$r5a$a$k(B. cf. liftWeight
1858: %>
1859: /flipWeight {
1860: /arg3 set /arg2 set /arg1 set
1861: [/cone /facet_i /ep /vp /v /v /ii] pushVariables
1862: [
1863: /cone arg1 def /facet_i arg2 def
1864: facet_i to_int32 /facet_i set
1865: /ep arg3 def
1866:
1867: ep to_univNum (1).. div /ep set
1868:
1869: % note: 2004.9.2
1870: cone (facetsv) getNode 2 get facet_i get /v set
1871: cone (facets) getNode 2 get facet_i get /f set
1872: /vp v 0 get def
1873: 1 1 v length 1 sub {
1874: /ii set
1875: vp v ii get add /vp set
1876: } for
1877: vp ep f mul sub /vp set
1878: vp nnormalize_vec /vp set
1879: /arg1 vp def
1880: ] pop
1881: popVariables
1882: arg1
1883: } def
1884:
1885: %<
1886: % Usages: cone1 cone2 isSameCone bool
1887: % cone1 cone2 $B$,Ey$7$$$+(B? facet $B$GHf$Y$k(B.
1888: % cone1, cone2 $B$O(B pointed cone $B$G$J$$$H$$$1$J$$(B.
1889: %>
1890: /isSameCone {
1891: /arg2 set /arg1 set
1892: [/cone1 /cone2 /facets1 /facets2 /ans] pushVariables
1893: [
1894: /cone1 arg1 def
1895: /cone2 arg2 def
1896: /facets1 cone1 (facets) getNode 2 get def
1897: /facets2 cone2 (facets) getNode 2 get def
1898: facets1 length facets2 length eq {
1899: facets1 facets2 sub isZero /ans set
1900: } {
1901: /ans 0 def
1902: } ifelse
1903: /arg1 ans def
1904: ] pop
1905: popVariables
1906: arg1
1907: } def
1908:
1909: %<
1910: % Usages: cone1 cone2 getCommonFacet list
1911: % cone1 $B$NCf$G(B cone2 $B$K4^$^$l$k(B facet $B$N%j%9%H(B
1912: % cone2 $B$NCf$G(B cone1 $B$K4^$^$l$k(B facet $B$N%j%9%H$r$b$I$9(B.
1913: % [1 [i] [j]] $B$"$k$H$-(B. [0 [ ] [ ]] $B$J$$$H$-(B.
1914: % cone1 $B$N(B facetsv[i] $B$,(B cone2 $B$K4^$^$l$k$+D4$Y$k(B.
1915: % cone2 $B$N(B facetsv[i] $B$,(B cone1 $B$K4^$^$l$k$+D4$Y$k(B.
1916: % cone1, cone2 $B$O(B pointed cone $B$G$J$$$H$$$1$J$$(B.
1917: %>
1918: /getCommonFacet {
1919: /arg2 set /arg1 set
1920: [/cone1 /cone2 /facets /ineq /ans1 /ans2 /i /tt] pushVariables
1921: [
1922: /cone1 arg1 def
1923: /cone2 arg2 def
1924:
1925: /facets cone1 (facetsv) getNode 2 get def
1926: /ineq cone2 (inequalities) getNode 2 get def
1927: /ans1 [
1928: 0 1 facets length 1 sub {
1929: /i set
1930: facets i get /tt set % facetsv[i] $B$r(B tt $B$X(B.
1931: ineq tt transpose mul isNonNegative {
1932: i
1933: } { } ifelse
1934: } for
1935: ] def
1936:
1937: /facets cone2 (facetsv) getNode 2 get def
1938: /ineq cone1 (inequalities) getNode 2 get def
1939: /ans2 [
1940: 0 1 facets length 1 sub {
1941: /i set
1942: facets i get /tt set % facetsv[i] $B$r(B tt $B$X(B.
1943: ineq tt transpose mul isNonNegative {
1944: i
1945: } { } ifelse
1946: } for
1947: ] def
1948: ans1 length 1 gt ans2 length 1 gt or {
1949: (getCommonFacet found more than 1 common facets.) error
1950: } { } ifelse
1951: % $B6&DL(B facet $B$,$"$l$P(B 1, $B$J$1$l$P(B 0.
1952: ans1 length 1 eq ans2 length 1 eq and {
1953: /tt 1 def
1954: } {
1955: /tt 0 def
1956: } ifelse
1957: /arg1 [tt ans1 ans2] def
1958: ] pop
1959: popVariables
1960: arg1
1961: } def
1962:
1963: %
1964: % -------------------------------------------------
1965: % test8 $B$O(B aux-cone.sm1 $B$X0\F0(B.
1966: % $B0J2<$$$h$$$h0lHL$N%W%m%0%i%`$N:n@.3+;O(B.
1967: % -------------------------------------------------
1968: %
1969:
1970: %<
1971: % Usages: setWeightBorder
1972: % cone.weightBorder (weight cone $B$N(B facet $B%Y%/%H%k$N=89g(B) $B$r@_Dj$9$k(B.
1973: % $B$"$HI{;:J*$H$7$F(B cone.w_cone_projectedWt (doPolymakeObj)
1974: % cone.w_ineq_projectedWt
1975: % cone.m $B<!85$N%Y%/%H%k(B.
1976: % cone.W, cone.Wt, cone.w_ineq $B$,$9$G$K7W;;$:$_$G$J$$$H$$$1$J$$(B.
1977: %>
1978: /setWeightBorder {
1979: [
1980: (Entering setWeightBorder ) message
1981: cone.w_ineq cone.Wt mul pruneZeroVector /cone.w_ineq_projectedWt set
1982: {
1983: cone.w_ineq_projectedWt length 0 eq {
1984: % weight $B$N6u4V$K(B border $B$,$J$$>l9g(B.
1985: /cone.weightBorder [ ] def
1986: exit
1987: } { } ifelse
1988: % weight $B$N6u4V$K(B border $B$,$"$k>l9g(B.
1989: cone.w_ineq_projectedWt getConeInfo /cone.w_cone_projectedWt set
1990: cone.w_cone_projectedWt 0 get 0 get to_int32 cone.m to_int32 eq {
1991: } {
1992: (setWeightBorder : internal error.) message
1993: } ifelse
1994: cone.w_cone_projectedWt 1 get (FACETS) getNode 2 get 0 get
1995: removeFirstFromPolymake /cone.weightBorder set
1996: exit
1997: } loop
1998: (cone.weightBorder=) message
1999: cone.weightBorder pmat
2000: ] pop
2001: } def
2002:
2003: %
2004: % -------------------------------------------------
2005: % $B%W%m%0%i%`$NN.$l(B.
2006: % Global: cone.fan cone $B$rG[Ns$H$7$F3JG<$9$k(B.
2007: %
2008: % ncone (next cone) $B$,?75,$KF@$i$l$?(B cone $B$G$"$k$H$9$k(B.
2009: % $B$3$N$H$-<!$NA`:n$r$9$k(B.
2010: % 0. ncone $B$,(B cone.fan $B$K$9$G$K$J$$$+D4$Y$k(B. $B$"$l$P(B, internal error.
2011: % 1. ncone markBorder ; ncone $B$NCf$N(B border $B>e$N(B facet $B$r(B mark
2012: % 2. cone.fan $B$NCf$N(B cone $B$H6&DL(B facet $B$,$J$$$+D4$Y(B (getCommonFacet),
2013: % $B$"$l$P$=$l$i$r(B mark $B$9$k(B.
2014: % global: cone.incidence $B$K(B $B6&DL(Bfacet $B$r;}$DAH$_$N>pJs$r2C$($k(B.
2015: % 3. ncone $B$r(B cone.fan $B$N:G8e$K2C$($k(B.
2016: % $B0J>e$NA`:n$r$^$H$a$?$b$N$,(B ncone updateFan
2017: %
2018: % getNextFlip $B$O(B cone.fan $B$NCf$+$i(B flip $B$7$F$J$$(B cone $B$H(B facet $B$NAH$rLa$9(B.
2019: % $B$J$1$l$P(B null $B$rLa$9(B. null $B$,La$l$P%W%m%0%i%`=*N;(B.
2020: %
2021: % getStargingCone $B$O7W;;$r=PH/$9$Y$-?75,$N(B cone $B$r7W;;$9$k(B. $BBg0hJQ?t(B cone.Lt, cone.W
2022: % $B$J$I$b$3$NCf$G@_Dj$9$k(B.
2023: % $BJQ?t%j%9%H(B, weight space $B$r@8@.$9$k4X?t(B, $BF~NOB?9`<0(B, weight $B$N8uJd(B $BEy$OBg0hJQ?t(B
2024: % $B$H$7$FF~NO$7$F$*$/(B.
2025: %
2026: % reduced gb $B$O(B $B4X?t(B input weight cone.gb reduced_G $B$G7W;;$9$k(B.
2027: %
2028: %
2029: % [ccone i] getNextCone ncone : flip $B$K$h$j<!$N(B cone $B$rF@$k(B.
2030: %
2031: % 1. clearGlobals ; $BF~NOBg0hJQ?t$N@_Dj(B.
2032: % 2. getStartingCone /ncone set
2033: % 3. { ncone updateFan
2034: % 4. getNextFlip /cone.nextflip set
2035: % 6. cone.nextflip isNull { exit } { } ifelse
2036: % 7. cone.nextflip getNextCone /ncone set
2037: % 8. } loop
2038: %
2039: %
2040: % -------------------------------------------------
2041: %
2042:
2043: %<
2044: % Usages: input weight cone.gb_Dh reduced_G
2045: % gb in h[1,1](D)
2046: %>
2047: /cone.gb_Dh {
2048: /arg2 set /arg1 set
2049: [/ff /ww /gg] pushVariables
2050: [
2051: /ff arg1 def
2052: /ww arg2 def
2053: [(AutoReduce) 1] system_variable
2054: [cone.vv ring_of_differential_operators
2055: [ww] weight_vector 0] define_ring
1.9 takayama 2056: [ff {toString .} map] ff getAttributeList setAttributeList
2057: groebner 0 get /gg set
1.1 takayama 2058: /cone.gb_Dh.g gg def
2059: /arg1 gg def
2060: ] pop
2061: popVariables
2062: arg1
2063: } def
2064:
2065: %<
2066: % Usages: cone.boundp
2067: %
2068: /cone.boundp {
2069: dup boundp 2 1 roll tag 0 eq not and
2070: } def
2071:
2072: %<
2073: % Usages: clearGlobals
2074: % cf. cone.boundp
2075: % polymake $B$r:FEY8F$V$?$a$K(B global $BJQ?t$r%/%j%"$9$k(B.
2076: % $B$^$@ESCf(B.
2077: %>
2078: /clearGlobals {
2079: /cone.W null def
2080: /cone.Wt null def
2081:
2082: /cone.cinit null def
2083: /cone.weightBorder null def
2084:
2085: } def
2086:
2087: %<
2088: % Usages: getStartingCone ncone
2089: % getStargingCone $B$O7W;;$r=PH/$9$Y$-?75,$N(B cone $B$r7W;;$9$k(B.
2090: % $B@_Dj$9$Y$-Bg0hJQ?t$O0J2<$r8+$h(B.
2091: %>
2092:
2093: /getStartingCone.test {
2094: %------------------Globals----------------------------------------
2095: % --------------- $BF~NO%G!<%?MQBg0hJQ?t$N@_Dj(B --------------------------
2096: %
2097: % cone.input : $BF~NOB?9`<07O(B
2098: /cone.input
2099: [(t1-x-y) (h*t2-x^2-y^2) (2*x*Dt2+h*Dt1+h*Dx) (2*y*Dt2+h*Dt1+h*Dy)]
2100: def
2101:
2102: % cone.vlist : $BA4JQ?t$N%j%9%H(B
2103: /cone.vlist [(t1) (t2) (x) (y) (Dt1) (Dt2) (Dx) (Dy) (h)] def
2104:
2105: % cone.vv : define_ring $B7A<0$NJQ?t%j%9%H(B.
2106: % t1,t2, x,y : t-space $B$N(B Grobner fan (local) $B$r5a$a$k(B.
2107: /cone.vv (t1,t2,x,y) def
2108:
2109: % cone.parametrizeWeightSpace : weight $B6u4V$r(B parametrize $B$9$k4X?t(B.
2110: % $BBg0hJQ?t(B cone.W , cone.Wpos $B$b$-$^$k(B.
2111: /cone.parametrizeWeightSpace {
2112: 4 2 parametrizeSmallFan
2113: } def
2114:
2115: % cone.w_start : weight$B6u4V$K$*$1$k(B weight $B$N=i4|CM(B.
2116: % $B$3$NCM$G(B max dim cone $B$,F@$i$l$J$$$H(B random weight $B$K$h$k(B $B%5!<%A$,;O$^$k(B.
2117: /cone.w_start
2118: [ 1 4 ]
2119: def
2120:
2121: % cone.gb : gb $B$r7W;;$9$k4X?t(B.
2122: /cone.gb {
2123: cone.gb_Dh
2124: } def
2125:
2126: %
2127: % ----------------- $B$*$o$j(B ---------------------------
2128: %
2129: } def % end of getStartingCone.test
2130:
2131: /getStartingCone {
2132: [/wv_start /w_start /reduced_G] pushVariables
2133: [
2134: % cone.n $B$O<+F0E*$K$-$a$i$l$k(B.
2135: % cone.n $B$O(B GB $B$r7W;;$9$k6u4V$N<!85(B.
2136: /cone.n cone.vlist length def
2137: %[1] cone.W, cone.Wpos $B$r5a$a$k(B. cone.m $B$O(B cone.W $B$h$j<+F0E*$K$-$^$k(B.
2138: % cone.m $B$O(B weight $B6u4V$N<+M3EY(B. cone.W $B$G<M1F$5$l$k@h$N<!85(B.
2139: /cone.W cone.boundp {
2140: (Skip cone.parametrizeWeightSpace. cf. clearGlobals) message
2141: } {
2142: cone.parametrizeWeightSpace
2143: } ifelse
2144: (parametrizing weight space: cone.W = ) messagen cone.W message
2145: /cone.Wt cone.W transpose def
2146: /cone.m cone.W length def
2147: % WeightBorder $B$N>r7oH=Dj(B facet $B$r@_Dj(B.
2148: /cone.weightBorder cone.boundp {
2149: (Skip setWeightBorder cf. clearGlobals) message
2150: } {
2151: setWeightBorder
2152: } ifelse
2153:
2154: %[2] weight vector wv_start $B$r@8@.$9$k(B.
2155: % wv_start $B$r@_Dj(B.
2156: cone.w_start tag 0 eq {
2157: % cone.w_start $B$,(B null $B$J$i(B random $B$K(B weight $B$r@_Dj(B.
2158: /cone.w_start cone.m cone_random_vec def
2159: } {
2160: cone.w_start length cone.m to_int32 eq {
2161: } {
2162: (Error: cone.w_start has wrong length.) error
2163: /cone.w_start cone.m cone_random_vec def
2164: } ifelse
2165: } ifelse
2166: /w_start cone.w_start cone.W mul def
2167:
2168: {
2169: cone.vlist w_start cone_wtowv /wv_start set
2170: (Trying a starting weight vector : ) messagen
2171: wv_start pmat
2172: %[3] reduced GB $B$N7W;;(B.
2173: cone.input wv_start cone.gb /reduced_G set
1.2 takayama 2174: (Reduced GB is obtained: ) message
2175: %reduced_G pmat
2176: /cone.cgb reduced_G def
2177: [cone.w_start w_start wv_start] /cone.cgb_weight set
1.1 takayama 2178:
2179: %[4] $B<M1F$7$F$+$i(B polytope $B$N%G!<%?$r7W;;(B.
2180: wv_start reduced_G coneEq /cone.g_ineq set
2181: cone.g_ineq cone.w_ineq join /cone.gw_ineq set
2182: cone.gw_ineq cone.Wt mul /cone.gw_ineq_projectedWt set % $B<M1F(B
2183: /cone.cinit cone.boundp {
2184: (Skipping cone.gw_ineq_projectedWt getConeInfo. cf. clearGlobals) message
2185: } {
2186: cone.gw_ineq_projectedWt getConeInfo /cone.cinit set
2187: } ifelse
2188:
2189: (cone.cinit is --- the first number is the dim of cone.) messagen
2190: cone.cinit 0 get pmat
2191: % Maximal dimensional cone $B$+$I$&$+$N8!::(B. $B8!::$K%Q%9$9$l$P(B loop $B$r(B exit
2192: % $B%Q%9$7$J$$>l9g(B w_start $B$r(B cone_random_vec $B$rMQ$$$FJQ99$9$k(B.
2193: cone.cinit 0 get 0 get to_int32 cone.m eq { exit }
2194: {
2195: (Failed to get the max dim cone. Updating the weight ...) messagen
1.2 takayama 2196: cone.m cone_random_vec /cone.w_start set
2197: /w_start cone.w_start cone.W mul def
1.1 takayama 2198: % cone.cinit $B$r:FEY7W;;$9$k$?$a$K(B clear $B$9$k(B.
2199: /cone.cinit null def
2200: } ifelse
2201: } loop
2202:
2203: (cone.m = ) messagen cone.m message
2204: (Suceeded to get the maximal dimensional startingCone.) message
2205:
2206: % Linearity subspace $B$N(B orth complement $B$X$N<M1F9TNs(B.
2207: % $BBg0hJQ?t(B cone.Lp, cone.Lpt $B$r@_Dj(B
2208: cone.cinit 0 get 1 get /cone.Lp set
2209: cone.Lp transpose /cone.Lpt set
2210: % Linearity subspace $B$N9TNs$r@_Dj(B.
2211: % $BBg0hJQ?t(B cone.L $B$r@_Dj(B
2212: cone.cinit 0 get 2 get /cone.L set
2213: % cone.d $B$O(B cone.W $B$*$h$S(B Linearity space $B$G3d$C$?8e(B, cone $B$r9M$($k$H$-$N<!85(B.
2214: % $BBg0hJQ?t(B cone.d $B$N@_Dj(B.
2215: /cone.d cone.Lp length def
2216:
2217: cone.m cone.d eq {
2218: (There is no linearity space) message
2219: } {
2220: (Dim of the linearity space is ) messagen cone.m cone.d sub message
2221: (cone.Lp = ) messagen cone.Lp pmat
2222: } ifelse
2223:
2224: %[5] cone.g_ineq * cone.Wt * cone.Lpt
2225: % cone.w_ineq * cone.Wt * cone.Lpt
2226: % $B$G@)Ls$r(B d $B<!85%Y%/%H%k$KJQ49(B.
2227: % W (R^m) $B6u4V$NITEy<0@)Ls$r(B L' (R^d) $B6u4V$X<M1F(B
2228: % cone.gw_ineq_projectedWtLpt
2229: % = cone.g_ineq*cone.Wt*cone.Lpt \/ cone.w_ineq*coneWt*cone.Lpt
2230:
2231: /cone.gw_ineq_projectedWtLpt
2232: cone.gw_ineq_projectedWt cone.Lpt mul
2233: def
2234:
2235: cone.m cone.d eq {
2236: /cone.cinit.d cone.cinit def
2237: } {
2238: % cone.m > cone.d $B$J$i$P(B, $B:FEY(B cone $B$N7W;;$,I,MW(B.
2239: % R^d $B$N(B cone $B$O(B cone.cinit.d $B$XF~$l$k(B.
2240: cone.gw_ineq_projectedWtLpt getConeInfo /cone.cinit.d set
2241: } ifelse
2242:
2243: cone.cinit.d 1 get newCone /cone.startingCone set
2244:
2245: (cone.startingCone is ) message
2246: cone.startingCone message
2247: ] pop
2248: popVariables
2249: cone.startingCone
2250: } def
2251:
2252: %
2253: % data/test9.sm1 $B$N(B test9 1-simplex X 2-simplex
2254: %
2255: % data/test10.sm1 1-simplex X 3-simplex
2256: % data/test11.sm1 SST, p.59
2257: %
2258: % $B$$$h$$$h(B, cone enumeration $B$N%W%m%0%i%`=q$-3+;O(B
2259: %
2260:
2261: %<
2262: % Usages: cone markBorder
2263: % cone->facets[i] $B$,(B weight space $B$N(B border $B$K$"$k$H$-(B
2264: % cone->flipped[i] = 2 $B$H$9$k(B.
2265: % $B$3$l$r(B cone $B$N$9$Y$F$N(B facet $B$KBP$7$F7W;;(B.
2266: %>
2267: /markBorder {
2268: /arg1 set
1.4 takayama 2269: [/cone /facets_t /flipped_t /kk /nextcid_t /nextfid_t] pushVariables
1.1 takayama 2270: [
2271: /cone arg1 def
2272: cone (facets) getNode 2 get /facets_t set
2273: cone (flipped) getNode 2 get /flipped_t set
1.4 takayama 2274: cone (nextcid) getNode 2 get /nextcid_t set
2275: cone (nextfid) getNode 2 get /nextfid_t set
1.1 takayama 2276: 0 1 flipped_t length 1 sub {
2277: /kk set
2278: flipped_t kk get (0).. eq {
2279: cone kk isOnWeightBorder {
2280: % Border $B$N>e$K$"$k$N$G(B flip $B:Q$N%^!<%/$r$D$1$k(B.
2281: flipped_t kk (2).. put
1.4 takayama 2282: % $B$H$J$j$N(B cone $B$N(B id (nextcid, nextfid) $B$O(B -2 $B$H$9$k(B.
2283: nextcid_t kk (-2).. put
2284: nextfid_t kk (-2).. put
1.1 takayama 2285: } { } ifelse
2286: } { } ifelse
2287: } for
2288: ] pop
2289: popVariables
2290: } def
2291:
2292: %<
2293: % Usages: ncone updateFan
2294: % $B%0%m!<%P%kJQ?t(B cone.fan $B$r99?7$9$k(B.
2295: %>
2296: %
2297: % updateFan $B$N(B debug $B$O(B data/test8 $B$G$H$j$"$($:$d$k(B.
2298: % test8 /ncone set $B$r<B9T$7$F$+$i(B ncone updateFan
2299:
2300: % global: cone.fan
2301: /cone.fan [ ] def
2302: % global: cone.incidence
2303: /cone.incidence [ ] def
1.2 takayama 2304: % global: cone.gblist gb's standing for each cones in cone.fan.
2305: /cone.gblist [ ] def
1.1 takayama 2306:
2307: /updateFan {
2308: /arg1 set
2309: [/ncone /kk /cfacet /ii /jj /tcone /flipped_t] pushVariables
2310: [
2311: /ncone arg1 def
2312: /cone.fan.n cone.fan length def
1.2 takayama 2313: % -1. cone.cgb ($BD>A0$K7W;;$5$l$?(B gb) $B$H(B cone.cgb_weight ($BD>A0$N7W;;$N(B weight)
2314: % $B$r(B cone.gblist $B$X3JG<$9$k(B.
2315: cone.gblist [ [cone.cgb cone.cgb_weight] newConeGB ] join /cone.gblist set
1.1 takayama 2316: % 0. ncone $B$,(B cone.fan $B$K$9$G$K$"$l$P%(%i!<(B
2317: 0 1 cone.fan.n 1 sub {
2318: /kk set
2319: ncone cone.fan kk get isSameCone {
2320: (Internal error updateFan: ncone is already in cone.fan) error
2321: } { } ifelse
2322: } for
2323:
2324: % 1. ncone $B$NCf$N(B border $B>e$N(B facet $B$r$9$Y$F(B mark.
2325: ncone markBorder
2326:
2327: % 2. ncone /\ cone.fan[kk] $B$,$"$k$+D4$Y$k(B. $B$"$l$P(B Mark $B$9$k(B. incidence graph $B$K2C$($k(B
2328: 0 1 cone.fan.n 1 sub {
2329: /kk set
2330: ncone cone.fan kk get getCommonFacet /cfacet set
2331: cfacet 0 get
2332: {
2333: % $B6&DL(B facet $B$,$"$k>l9g(B. [[cone$BHV9f(B face$BHV9f(B] [cone$BHV9f(B face$BHV9f(B]] $B$N7A<0$G3JG<(B.
2334: /ii cfacet 1 get 0 get def
2335: /jj cfacet 2 get 0 get def
2336: cone.incidence [ [[cone.fan.n ii] [kk jj]] ] join /cone.incidence set
2337: % flipped $B$r(B mark $B$9$k(B.
2338: ncone ii markFlipped
2339: cone.fan kk get /tcone set
2340: tcone jj markFlipped
1.4 takayama 2341: % nextcid, nextfid $B$r@_Dj$9$k(B.
2342: ncone ii [kk jj] markNext
2343: tcone jj [cone.fan.n ii] markNext
1.1 takayama 2344: } { } ifelse
2345: } for
2346: % 3. ncone $B$r2C$($k(B.
2347: cone.fan [ncone] join /cone.fan set
2348: ] pop
2349: popVariables
2350: } def
2351:
2352: %<
1.9 takayama 2353: % usages: getNextFlip [cone, k, cid]
1.1 takayama 2354: % cone.fan $B$r8!:w$7$F(B $B$^$@(B flip $B$7$F$J$$(B cone $B$H(B facet $B$NAH$rLa$9(B.
2355: % $B$b$&$J$$$H$-$K$O(B null $B$rLa$9(B.
1.9 takayama 2356: % cid $B$O(B cone $B$,(B cone.fan $B$N(B $B2?HVL\$G$"$k$+$N(B index. cone.gblist $B$N8!:wEy$K(B
2357: % $BMQ$$$k(B.
1.1 takayama 2358: %>
2359: /getNextFlip {
1.9 takayama 2360: [/tcone /ans /ii /cid] pushVariables
1.1 takayama 2361: [
1.9 takayama 2362: /ans null def /cid -1 def
1.1 takayama 2363: 0 1 cone.fan length 1 sub {
2364: /ii set
2365: cone.fan ii get /tcone set
1.9 takayama 2366: /cid ii def
1.1 takayama 2367: tcone getNextFacet /ans set
2368: ans tag 0 eq { } { exit } ifelse
2369: } for
2370: ans tag 0 eq { /arg1 null def }
1.9 takayama 2371: { /arg1 [tcone ans cid] def } ifelse
1.1 takayama 2372: ] pop
2373: popVariables
2374: arg1
2375: } def
2376:
2377: % global variable : cone.epsilon , cone.epsilon.limit
2378: % flip $B$N;~$N(B epsilon
2379: /cone.epsilon (1).. (10).. div def
2380: /cone.epsilon.limit (1).. (100).. div def
1.9 takayama 2381: % cone.epsilon.limit $B$rIi$K$9$l$PDd;_$7$J$$(B.
1.1 takayama 2382:
2383: %<
2384: % Usages: result_getNextFlip getNextCone ncone
2385: % flip $B$7$F?7$7$$(B ncone $B$rF@$k(B.
2386: %>
2387: /getNextCone {
2388: /arg1 set
2389: [/ncone /ccone /kk /w /next_weight_w_wv] pushVariables
2390: [
2391: /ccone arg1 def
2392: /ncone null def
2393: /kk ccone 1 get def
2394: ccone 0 get /ccone set
2395: {
2396: ccone tag 0 eq { exit } { } ifelse
2397:
2398: % ccone $B$N(B kk $BHVL\$N(B facet $B$K$D$$$F(B flip $B$9$k(B.
2399: ccone kk cone.epsilon flipWeight /w set
2400: (Trying new weight is ) messagen w message
2401: w liftWeight /next_weight_w_wv set
2402: (Trying new weight [w,wv] is ) messagen next_weight_w_wv message
2403:
2404: cone.input next_weight_w_wv 1 get cone.gb /cone.cgb set
1.2 takayama 2405: [w] next_weight_w_wv join /cone.cgb_weight set
1.1 takayama 2406: next_weight_w_wv 1 get cone.cgb coneEq /cone.g_ineq set
2407: cone.g_ineq cone.w_ineq join cone.Wt mul cone.Lpt mul
2408: pruneZeroVector /cone.gw_ineq_projectedWtLpt set
2409:
2410: (cone.gw_ineq_projectedWtLpt is obtained.) message
2411:
2412: cone.gw_ineq_projectedWtLpt getConeInfo /cone.nextConeInfo set
2413: % $B<!85$rD4$Y$k(B. $B$@$a$J$i(B retry
2414: cone.nextConeInfo 0 get 0 get to_int32 cone.d eq {
2415: cone.nextConeInfo 1 get newCone /ncone set
2416: ccone ncone getCommonFacet 0 get {
2417: (Flip succeeded.) message
2418: exit
2419: } { } ifelse
2420: } { } ifelse
2421: % common face $B$,$J$1$l$P(B $B$d$O$j(B epsilon $B$r>.$5$/(B.
2422: cone.nextConeInfo 0 get 0 get to_int32 cone.d eq {
2423: (ccone and ncone do not have a common facet.) message
2424: } {
2425: (ncone is not maximal dimensional. ) message
2426: } ifelse
2427: (Decreasing epsilon to ) messagen
2428: cone.epsilon (1).. (2).. div mul /cone.epsilon set
2429: cone.epsilon cone.epsilon.limit sub numerator (0).. lt {
2430: (Too small cone.epsilon ) error
2431: } { } ifelse
2432: cone.epsilon message
2433: } loop
2434: /arg1 ncone def
2435: ] pop
2436: popVariables
2437: arg1
2438: } def
2439:
2440: %<
2441: % Usages: set globals and getGrobnerFan
2442: % cf. clearGlobals
2443: % getStartingCone $B$9$k$H(B weightSpace $B$H$+$N7W;;$,$G$-$k(B. isOnWeightBorder $B$,(B
2444: % $B7h$a$i$l$k(B.
2445: %>
2446: % $B$H$j$"$($:(B (data/test8.sm1) run $B$7$F$+$i(B getGrobnerFan
2447: /getGrobnerFan {
2448: getStartingCone /cone.ncone set
2449: {
2450: cone.ncone updateFan
2451: ( ) message
2452: (----------------------------------------------------------) message
2453: (getGrobnerFan #cone.fan=) messagen cone.fan length message
2454: cone.ncone /cone.ccone set
2455: getNextFlip /cone.nextflip set
2456: cone.nextflip tag 0 eq { exit } { } ifelse
2457: cone.nextflip getNextCone /cone.ncone set
2458: } loop
1.2 takayama 2459: (Construction is completed. See cone.fan, cone.incidence and cone.gblist.)
2460: message
2461: } def
2462:
2463: %<
2464: % Usages: vlist generateD1_1
2465: % -1,1 weight $B$r@8@.$9$k(B.
2466: % vlist $B$O(B (t,x,y) $B$+(B [(t) (x) (y)]
2467: %
2468: %>
2469: /generateD1_1 {
2470: /arg1 set
2471: [/vlist /rr /rr /ii /vv] pushVariables
2472: [
2473: /vlist arg1 def
2474: vlist isString {
2475: [vlist to_records pop] /vlist set
2476: } { } ifelse
2477: [
2478: 0 1 vlist length 1 sub {
2479: /ii set
2480: vlist ii get /vv set
2481: vv -1
2482: [@@@.Dsymbol vv] cat 1
2483: } for
2484: ] /rr set
2485: /arg1 rr def
2486: ] pop
2487: popVariables
2488: arg1
2489: } def
2490:
2491: /listNodes {
2492: /arg1 set
2493: [/in-listNodes /ob /rr /rr /ii] pushVariables
2494: [
2495: /ob arg1 def
2496: /rr [ ] def
2497: {
2498: ob isClass {
2499: ob (array) dc /ob set
2500: } { exit } ifelse
2501: rr [ob 0 get] join /rr set
2502: ob 2 get /ob set
2503: 0 1 ob length 1 sub {
2504: /ii set
2505: rr ob ii get listNodes join /rr set
2506: } for
2507: exit
2508: } loop
2509: /arg1 rr def
2510: ] pop
2511: popVariables
2512: arg1
2513: } def
2514: [(listNodes)
2515: [(ob listNodes)
2516: (cf. getNode)
2517: (Example:)
2518: ( /dog [(dog) [[(legs) 4] ] [ ]] [(class) (tree)] dc def)
2519: ( /man [(man) [[(legs) 2] ] [ ]] [(class) (tree)] dc def)
2520: ( /ma [(mammal) [ ] [man dog]] [(class) (tree)] dc def)
2521: ( ma listNodes )
2522: ]] putUsages
2523:
2524: %<
2525: % Usages: obj printTree
2526: %>
2527: /printTree {
2528: /arg1 set
2529: [/ob /rr /rr /ii /keys /tt] pushVariables
2530: [
2531: /ob arg1 def
2532: /rr [ ] def
2533: /keys ob listNodes def
2534: keys 0 get /tt set
2535: keys rest /keys set
2536: keys { ob 2 1 roll getNode } map /rr set
2537: (begin ) messagen tt messagen
2538: ( ---------------------------------------) message
2539: 0 1 rr length 1 sub {
2540: /ii set
2541: keys ii get messagen (=) message
2542: rr ii get 2 get pmat
2543: } for
2544: (--------------------------------------- end ) messagen
2545: tt message
2546: /arg1 rr def
2547: ] pop
2548: popVariables
2549: arg1
2550: } def
2551:
2552: %<
2553: % Usages $B$O(B (inputForm) usages $B$r$_$h(B.
2554: %>
2555: /inputForm {
2556: /arg1 set
2557: [/ob /rr /i ] pushVariables
2558: [
2559: /ob arg1 def
2560: /rr [ ] def
2561: {
2562: ob isArray {
2563: rr [ ([) ] join /rr set
2564: 0 1 ob length 1 sub {
2565: /i set
2566: i ob length 1 sub lt {
2567: rr [ob i get inputForm $ , $] join /rr set
2568: } {
2569: rr [ob i get inputForm] join /rr set
2570: } ifelse
2571: } for
2572: rr [ (]) ] join cat /rr set
2573: exit
2574: } { } ifelse
2575: ob isClass {
2576: ob etag 263 eq { % tree
2577: /rr ob inputForm.tree def exit
2578: } { /rr [( $ this etag is not implemented $ )] cat def exit } ifelse
2579: } { } ifelse
2580: ob isUniversalNumber {
2581: [$($ ob toString $)..$] cat /rr set
2582: exit
2583: } { } ifelse
2584: ob isPolynomial {
2585: [$($ ob toString $).$] cat /rr set
2586: exit
2587: } { } ifelse
2588: ob isRational {
2589: [$ $ ob (numerator) dc inputForm $ $
2590: ob (denominator) dc inputForm $ div $ ] cat /rr set
2591: exit
2592: } { } ifelse
2593: ob isString {
2594: [$($ ob $)$ ] cat /rr set
2595: exit
2596: } { } ifelse
2597: ob toString /rr set
2598: exit
2599: } loop
2600: rr /arg1 set
2601: ] pop
2602: popVariables
2603: arg1
2604: } def
2605: [(inputForm)
2606: [(obj inputForm str)
2607: ]] putUsages
2608: % should be moved to dr.sm1
2609:
2610: /inputForm.tree {
2611: /arg1 set
2612: [/ob /key /rr /rr /ii] pushVariables
2613: [
2614: /ob arg1 def
2615: /rr [ ] def
2616: {
2617: ob (array) dc /ob set
2618: /rr [ $[$ ob 0 get inputForm $ , $
2619: ob 1 get inputForm $ , $
2620: ] def
2621: rr [ob 2 get inputForm ] join /rr set
2622: rr [$ ] $] join /rr set
2623: rr [ $ [(class) (tree)] dc $ ] join /rr set
2624: rr cat /rr set
2625: exit
2626: } loop
2627: /arg1 rr def
2628: ] pop
2629: popVariables
2630: arg1
2631: } def
2632:
2633: %<
2634: % Usages: str inputForm.value str
2635: %>
2636: /inputForm.value {
2637: /arg1 set
2638: [/key /val /valstr /rr] pushVariables
2639: [
2640: arg1 /key set
2641: key isString { } {(inputForm.value: argument must be a string) error } ifelse
2642: key boundp {
2643: [(parse) key] extension pop
2644: /val set
2645: val inputForm /valstr set
2646: [( ) valstr ( /) key ( set )] cat /rr set
2647: } {
2648: /valstr [] cat /rr set
2649: } ifelse
2650: rr /arg1 set
2651: ] pop
2652: popVariables
2653: arg1
2654: } def
2655:
2656: % global: cone.withGblist
2657: /cone.withGblist 0 def
2658: %<
2659: % Usages: saveGrobnerFan str
2660: % GrobnerFan $B$N%G!<%?$r(B inputForm $B$KJQ99$7$FJ8;zNs$KJQ$($k(B.
2661: % $B$3$N%G!<%?$r(B parse $B$9$k$H(B GrobnerFan $B$rF@$k$3$H$,2DG=(B.
2662: % BUG: $BB?9`<0$NB0$9$k4D$N%G!<%?$NJ]B8$O$^$@$7$F$J$$(B.
2663: %>
2664: /saveGrobnerFan {
2665: [/rr] pushVariables
2666: [
2667: (cone.withGblist=) messagen cone.withGblist message
2668: [
2669: % $B%f!<%6$N@_Dj$9$k%Q%i%a!<%?(B. cone.gb, cone.parametrizeWeightSpace $BEy$N4X?t$b$"$j(B.
2670: (cone.comment)
2671: (cone.type) (cone.local) (cone.h0)
2672: (cone.vlist) (cone.vv)
2673: (cone.input)
2674:
2675: % $B%W%m%0%i%`Cf$GMxMQ$9$k(B, $BBg;v$JBg0hJQ?t(B. weight vector $B$N<M1F9TNs$,=EMW(B.
2676: (cone.n) (cone.m) (cone.d)
2677: (cone.W) (cone.Wpos) (cone.Wt)
2678: (cone.L) (cone.Lp) (cone.Lpt)
2679: (cone.weightBorder)
2680: (cone.w_ineq)
2681: (cone.w_ineq_projectedWt)
2682: (cone.epsilon)
2683:
2684: % $B7k2L$NMWLs(B.
2685: (cone.fan)
2686: cone.withGblist { (cone.gblist) } { } ifelse
2687: (cone.incidence)
2688:
2689: ] { inputForm.value nl } map /rr set
1.3 takayama 2690: rr cat /rr set
2691: % ring $B$r(B save $B$7$F$J$$$N$GEv:B$NBP=h(B.
2692: [ ([) cone.vv inputForm ( ring_of_differential_operators 0 ] define_ring )
2693: nl nl rr] cat /arg1 set
1.2 takayama 2694: ] pop
2695: popVariables
2696: arg1
2697: } def
2698:
2699: /printGrobnerFan.1 {
2700: /arg1 set
2701: [/key /rr] pushVariables
2702: [
2703: /key arg1 def
2704: key boundp {
2705: [(parse) key] extension pop /rr set
2706: rr isArray {
2707: key messagen ( = ) message rr pmat
2708: } {
2709: key messagen ( = ) messagen rr message
2710: } ifelse
2711: }{
2712: key messagen ( = ) message
2713: } ifelse
2714: ] pop
2715: popVariables
2716: } def
2717:
2718: /printGrobnerFan {
2719: [/i] pushVariables
2720: [
2721: (========== Grobner Fan ====================) message
2722: [
2723: (cone.comment)
2724: (cone.vlist) (cone.vv)
2725: (cone.input)
2726: (cone.type) (cone.local) (cone.h0)
2727: (cone.n) (cone.m) (cone.d)
2728: (cone.W) (cone.Wpos) (cone.Wt)
2729: (cone.L) (cone.Lp) (cone.Lpt)
2730: (cone.weightBorder)
2731: (cone.incidence)
2732: ] { printGrobnerFan.1 } map
2733: ( ) message
2734: 0 1 cone.fan length 1 sub {
2735: /ii set
2736: ii messagen ( : ) messagen
2737: cone.fan ii get printTree
2738: } for
2739: cone.withGblist {
2740: 0 1 cone.gblist length 1 sub {
2741: /ii set
2742: ii messagen ( : ) messagen
2743: cone.gblist ii get printTree
2744: } for
2745: } { } ifelse
2746:
2747:
2748: (=========================================) message
2749: (cone.withGblist = ) messagen cone.withGblist message
2750: ( ) message
2751: ] pop
2752: popVariables
2753: } def
2754:
2755: %<
2756: % Usages: m uniq
2757: % Remove duplicated lines.
2758: %>
2759: /uniq {
2760: /arg1 set
2761: [/mm /prev /i /rr] pushVariables
2762: [
2763: /mm arg1 def
2764: {
2765: mm length 0 eq { [ ] /rr set exit } { } ifelse
2766: /prev mm 0 get def
2767: [
2768: prev
2769: 1 1 mm length 1 sub {
2770: /i set
2771: mm i get prev sub isZero { }
2772: { /prev mm i get def prev } ifelse
2773: } for
2774: ] /rr set
2775: exit
2776: } loop
2777: rr /arg1 set
2778: ] pop
2779: popVariables
2780: arg1
2781: } def
1.3 takayama 2782:
2783: %<
2784: % Usages: [vlist vw_vector] getGrRing [vlist vGlobal sublist]
2785: % example: [(x,y,z) [(x) -1 (Dx) 1 (y) 1 (Dy) 2]] getGrRing
2786: % [(x,y,z,y') [(x)] [[(Dy) (y')]]]
2787: % h[0,1](D_0) $B@lMQ$N(B getGrRing.
2788: % u_i + v_i > 0 $B$J$i(B Dx_i ==> x_i' ($B2D49$JJQ?t(B). sublist $B$X(B.
2789: % u_i < 0 $B$J$i(B x_i $B$O(B vGlobal $B$X(B.
2790: % ii [vlist vGlobal sublist] toGrRing /ii set
2791: % [ii jj vlist [(partialEcartGlobalVarX) vGlobal]] ecart.isSameIdeal $B$H;H$&(B.
2792: %>
2793: /getGrRing {
2794: /arg1 set
2795: [/vlist /vw_vector /ans /vGlobal /sublist /newvlist
2796: /dlist /tt /i /u /v /k
2797: ] pushVariables
2798: [
2799: /vlist arg1 0 get def
2800: /vw_vector arg1 1 get def
2801:
2802: vlist isString { [vlist to_records pop] /vlist set } { } ifelse
2803: vlist { toString } map /vlist set
2804: % dlist $B$O(B [(Dx) (Dy) (Dz)] $B$N%j%9%H(B.
2805: vlist { /tt set [@@@.Dsymbol tt] cat } map /dlist set
2806:
2807: /newvlist [ ] def /sublist [ ] def /vGlobal [ ] def
2808: % $B2D49$J?7$7$$JQ?t$r(B newvlist $B$X(B. $BCV49I=$r(B sublist $B$X(B.
2809: 0 1 vlist length 1 sub {
2810: /i set
2811: % (u,v) $B$O(B (x_i, Dx_i) $B$KBP$9$k(B weight vector
2812: /u vlist i get , vw_vector getGrRing.find def
2813: u -1 gt {
2814: vw_vector , u 1 add , get /u set
2815: } { /u 0 def } ifelse
2816:
2817: /v dlist i get , vw_vector getGrRing.find def
2818: v -1 gt {
2819: vw_vector , v 1 add , get /v set
2820: } { /v 0 def } ifelse
2821: u to_int32 /u set , v to_int32 /v set
2822:
2823: u v add , 0 gt {
2824: newvlist [vlist i get] join /newvlist set
2825: } { } ifelse
2826: u 0 lt {
2827: vGlobal [vlist i get] join /vGlobal set
2828: } { } ifelse
2829: } for
2830:
2831: newvlist { /tt set [ [@@@.Dsymbol tt] cat [tt (')] cat ] } map
2832: /sublist set
2833:
2834: /ans [ vlist , newvlist { /tt set [tt (')] cat } map , join from_records
2835: vGlobal sublist] def
2836: /arg1 ans def
2837: ] pop
2838: popVariables
2839: arg1
2840: } def
2841:
2842: %<
2843: % Usages: a uset getGrRing.find index
2844: %>
2845: /getGrRing.find {
2846: /arg2 set /arg1 set
2847: [/a /uset /ans /i] pushVariables
2848: [
2849: /a arg1 def /uset arg2 def
2850: /ans -1 def
2851: { /ans -1 def
2852: 0 1 , uset length 1 sub {
2853: /i set
2854: a tag , uset i get tag eq {
2855: a , uset i get eq {
2856: /ans i def exit
2857: } { } ifelse
2858: } { } ifelse
2859: } for
2860: exit
2861: } loop
2862: /arg1 ans def
2863: ] pop
2864: popVariables
2865: arg1
2866: } def
2867:
2868: %<
2869: % Usages: g1 g2 isSameGrRing bool
2870: % g1, g2 $B$O(B getGrRing $B$NLa$jCM(B.
2871: %>
2872: /isSameGrRing {
2873: /arg2 set /arg1 set
2874: [/g1 /g2 /ans] pushVariables
2875: [
2876: /g1 arg1 def /g2 arg2 def
2877: {
2878: /ans 1 def
2879: g1 0 get , g2 0 get eq { } { /ans 0 def exit } ifelse
2880: exit
2881: g1 1 get , g2 1 get eq { } { /ans 0 def exit } ifelse
2882: } loop
2883: /arg1 ans def
2884: ] pop
2885: popVariables
2886: arg1
2887: } def
2888:
2889: %<
2890: % Usages: [[ii i_vw_vector] [jj j_vw_vector] vlist] isSameInGrRing_h
1.4 takayama 2891: % It computes gb.
1.3 takayama 2892: %>
2893: /isSameInGrRing_h {
2894: /arg1 set
2895: [/ii /i_vw_vector /jj /j_vw_vector /vlist
2896: /i_gr /j_gr /rrule /ans] pushVariables
2897: [
2898: /ii arg1 [0 0] get def
2899: /i_vw_vector arg1 [0 1] get def
2900: /jj arg1 [1 0] get def
2901: /j_vw_vector arg1 [1 1] get def
2902: /vlist arg1 2 get def
2903: {
2904: [vlist i_vw_vector] getGrRing /i_gr set
2905: [vlist j_vw_vector] getGrRing /j_gr set
2906: i_gr j_gr isSameGrRing { } { /ans [0 [i_gr j_gr]] def exit} ifelse
2907:
2908: % bug: in case of module
2909: [i_gr 0 get , ring_of_differential_operators 0] define_ring
2910:
2911: % H $B$r(B 1 $B$K(B.
2912: /rrule [ [@@@.Hsymbol . (1).] ] def
2913:
2914: i_gr 2 get length 0 eq {
2915: } {
2916: rrule i_gr 2 get { { . } map } map join /rrule set
2917: } ifelse
2918: ii { toString . rrule replace toString } map /ii set
2919: jj { toString . rrule replace toString } map /jj set
2920:
2921: [ii jj i_gr 0 get , i_gr 1 get] ecartd.isSameIdeal_h /ans set
2922: [ans [i_gr] rrule ecartd.isSameIdeal_h.failed] /ans set
2923:
2924: exit
2925: } loop
2926: /arg1 ans def
2927: ] pop
2928: popVariables
2929: arg1
2930: } def
2931:
2932: /test1.isSameInGrRing_h {
2933: [(parse) (data/test8-data.sm1) pushfile] extension
2934:
2935: cone.gblist 0 get (initial) getNode 2 get /ii set
2936: cone.gblist 0 get (weight) getNode [2 0 2] get /iiw set
2937:
2938: cone.gblist 1 get (initial) getNode 2 get /jj set
2939: cone.gblist 1 get (weight) getNode [2 0 2] get /jjw set
2940:
2941: (Doing [ [ii iiw] [jj jjw] cone.vv ] isSameInGrRing_h /ff set) message
2942: [ [ii iiw] [jj jjw] cone.vv ] isSameInGrRing_h /ff set
2943:
2944: ff pmat
2945:
2946: } def
2947:
2948:
2949: %<
1.4 takayama 2950: % Usages: i j isSameCone_h.0 [bool, ...]
2951: % $B%F%9%HJ}K!(B. (data/test8.sm1) run (data/test8-data.sm1) run 0 1 isSameCone_h.0
2952: % gb $B$r:FEY7W;;$9$k(B stand alone $BHG(B. gr(Local ring) $B$GHf3S(B.
1.3 takayama 2953: %>
1.4 takayama 2954: /isSameCone_h.0 {
1.3 takayama 2955: /arg2 set /arg1 set
2956: [/i /j /ans /ii /iiw /jj /jjw] pushVariables
2957: [
2958: /i arg1 def /j arg2 def
1.4 takayama 2959: i to_int32 /i set , j to_int32 /j set
1.3 takayama 2960: cone.debug { (Comparing ) messagen [i j] message } { } ifelse
2961:
2962: cone.gblist i get (initial) getNode 2 get /ii set
2963: cone.gblist i get (weight) getNode [2 0 2] get /iiw set
2964:
2965: cone.gblist j get (initial) getNode 2 get /jj set
2966: cone.gblist j get (weight) getNode [2 0 2] get /jjw set
2967:
2968: [ [ii iiw] [jj jjw] cone.vv ] isSameInGrRing_h /ans set
2969:
2970: ans /arg1 set
2971: ] pop
2972: popVariables
2973: arg1
2974: } def
2975:
1.4 takayama 2976: %<
2977: % Usages: [ii vv i_vw_vector] getGbInGrRing_h [ii_gr i_gr]
2978: % Get Grobner Basis of ii in the graded ring.
2979: % The graded ring is obtained automatically from vv and i_vw_vector.
2980: % ii_gr is the Grobner basis. i_gr is the output of getGrRing.
2981: % cf. isSameInGrRing_h, ecart.isSameIdeal_h with [(noRecomputation) 1]
2982: %>
2983: /getGbInGrRing_h {
2984: /arg1 set
2985: [/ii /i_vw_vector /vlist /rng /vv /vvGlobal /wv /iigg
2986: /i_gr /rrule /ans] pushVariables
2987: [
2988: /ii arg1 0 get def
2989: /vlist arg1 1 get def
2990: /i_vw_vector arg1 2 get def
2991: [vlist i_vw_vector] getGrRing /i_gr set
2992:
2993: % bug: in case of module
2994: [i_gr 0 get , ring_of_differential_operators 0] define_ring
2995:
2996: % H $B$r(B 1 $B$K(B.
2997: /rrule [ [@@@.Hsymbol . (1).] ] def
2998:
2999: i_gr 2 get length 0 eq {
3000: } {
3001: rrule i_gr 2 get { { . } map } map join /rrule set
3002: } ifelse
3003: /vvGlobal i_gr 1 get def
3004: /vv i_gr 0 get def
3005:
3006: ii { toString . rrule replace toString } map /ii set
3007:
3008: [vv vvGlobal] ecart.stdBlockOrder /wv set
3009: vvGlobal length 0 eq {
3010: /rng [vv wv ] def
3011: }{
3012: /rng [vv wv [(partialEcartGlobalVarX) vvGlobal]] def
3013: } ifelse
3014: /save-cone.autoHomogenize ecart.autoHomogenize def
3015: /ecart.autoHomogenize 0 def
3016: [ii] rng join ecartd.gb /iigg set
3017: save-cone.autoHomogenize /ecart.autoHomogenize set
3018: /ans [iigg 0 get i_gr] def
3019: /arg1 ans def
3020: ] pop
3021: popVariables
3022: arg1
3023: } def
3024:
3025: /test1.getGbInGrRing_h {
3026: [(parse) (data/test8-data.sm1) pushfile] extension
3027:
3028: cone.gblist 0 get (initial) getNode 2 get /ii set
3029: cone.gblist 0 get (weight) getNode [2 0 2] get /iiw set
3030: [ii cone.vv iiw] getGbInGrRing_h /ff1 set
3031:
3032: cone.gblist 1 get (initial) getNode 2 get /jj set
3033: cone.gblist 1 get (weight) getNode [2 0 2] get /jjw set
3034: [jj cone.vv jjw] getGbInGrRing_h /ff2 set
3035:
3036: (ff1 and ff2) message
3037:
3038: } def
3039:
3040:
3041: %<
3042: % setGrGblist
3043: % cone.grGblist $B$r@_Dj$9$k(B.
3044: %>
3045: /setGrGblist {
3046: [/ii /ww /gg] pushVariables
3047: [
3048: cone.gblist {
3049: /gg set
3050: gg (initial) getNode 2 get /ii set
3051: gg (weight) getNode [2 0 2] get /ww set
3052: [ii cone.vv ww] getGbInGrRing_h
3053: } map /cone.grGblist set
3054: ] pop
3055: popVariables
3056: } def
3057:
3058: %<
3059: % Usages: i j isSameCone_h.2 [bool, ...]
3060: % gb $B$r:FEY7W;;$7$J$$(B.
3061: %>
3062: /isSameCone_h.2 {
3063: /arg2 set /arg1 set
3064: [/i /j /ans /ii /iiw /jj /jjw] pushVariables
3065: [
3066: /i arg1 def /j arg2 def
3067: i to_int32 /i set , j to_int32 /j set
3068: (cone.grGblist) boundp { } { setGrGblist } ifelse
3069: cone.debug { (Comparing ) messagen [i j] message } { } ifelse
3070:
3071: cone.grGblist i get /ii set
3072: cone.grGblist j get /jj set
3073:
3074: ii 1 get , jj 1 get isSameGrRing { }
3075: { /ans [0 [ii 1 get jj 1 get]] def exit} ifelse
3076:
3077: [ii 0 get , jj 0 get cone.vv [[(noRecomputation) 1]] ]
3078: ecartd.isSameIdeal_h /ans set
3079: [ans [ii 1 get] ii 1 get , ecartd.isSameIdeal_h.failed] /ans set
3080:
3081: ans /arg1 set
3082: ] pop
3083: popVariables
3084: arg1
3085: } def
3086:
3087: %<
3088: % test1.isSameCone_h.2 $B$O(B cone.grGblist $B$K(B initial $B$N(B gb $B$r(B graded ring
3089: % $B$G$^$:7W;;$7(B, $B$=$l$+$i(B ideal $B$NHf3S$r$*$3$J$&(B. isSameCone_h.1 $B$KHf$Y$F(B
3090: % gb $B$N:FEY$N7W;;$,$J$$$N$G7P:QE*(B.
3091: %>
3092: /test1.isSameCone_h.2 {
3093: /cone.loaded boundp { }
3094: {
3095: [(parse) (cohom.sm1) pushfile] extension
3096: [(parse) (dhecart.sm1) pushfile] extension
3097: /cone.loaded 1 def
3098: } ifelse
3099: %[(parse) (cone.sm1) pushfile] extension
3100: [(parse) (data/test8-data.sm1) pushfile] extension
3101: setGrGblist
3102: (cone.grGblist is set.) message
3103: 0 1 isSameCone_h.2 pmat
3104: } def
3105:
3106: %<
3107: % dhcone $B$O(B DeHomogenized Cone $B$NN,(B. H->1 $B$H$7$F(B cone $B$r(B merge $B$7$F$$$/4X?t(B
3108: % $B$dBg0hJQ?t$K;H$&(B.
3109: % cone.gblist, cone.fan $B$,@5$7$/@_Dj$5$l$F$$$k$3$H(B.
3110: % (setGrGblist $B$r<B9T:Q$G$"$k$3$H(B. $B<+F0<B9T$5$l$k$,(B... )
3111: %
3112: %>
3113:
3114: /isSameCone_h { isSameCone_h.2 } def
3115:
3116: %<
3117: % Usages: genDhcone.init
3118: % dhcone.checked (dehomogenized $B:Q$N(B cone$BHV9f(B), dhcone.unchecked $B$N=i4|2=(B.
3119: %>
3120: /genDhcone.init {
3121: /dhcone.checked [ ] def
3122: /dhcone.unchecked [
3123: 0 1 cone.fan length 1 sub {
3124: to_univNum
3125: } for
3126: ] def
3127: } def
3128:
3129: %<
3130: % Usages: k genDhcone dhcone
3131: % cone.fan[k] $B$r=PH/E@$H$7$F(B cone $B$r(B dehomogenize $B$9$k(B (merge $B$9$k(B).
3132: %
3133: % $B%F%9%H(B1. (data/test14.sm1) run (data/test14-data.sm1) run
3134: % genDhcone.init
3135: % 0 genDhcone /ff set
3136: %>
3137:
3138: /genDhcone {
3139: /arg1 set
3140: [/k /facets /merged /nextcid /nextfid /coneid
3141: /newfacets /newmerged /newnextcid /newnextfid /newconeid /vv
3142: /i /j /p /q /rr /cones /differentC
3143: ] pushVariables
3144: [
3145: /k arg1 def
3146: /facets [ ] def /merged [ ] def /nextcid [ ] def
3147: /nextfid [ ] def /coneid [ ] def
3148: /cones [ ] def
3149: /differentC [ ] def
3150:
3151: k to_univNum /k set
3152:
3153: {
3154: % Step1. cone.fan[k] $B$r(B $B2C$($k(B. new... $B$X=i4|%G!<%?$r=q$-9~$`(B.
3155: cone.debug {(Step 1. Adding ) messagen k messagen (-th cone.) message} { } ifelse
3156: cones [k to_univNum] join /cones set
3157: cone.fan k get , (facets) getNode 2 get /vv set
3158: /newfacets [ ] vv join def
3159:
3160: cone.fan k get , (nextcid) getNode 2 get /vv set
3161: /newnextcid [ ] vv join def
3162:
3163: cone.fan k get , (nextfid) getNode 2 get /vv set
3164: /newnextfid [ ] vv join def
3165:
3166: % newmerged $B$O$^$:(B 0 $B$G$&$a$k(B. 0 : $B$^$@D4$Y$F$J$$(B.
3167: % 1 : merged $B$G>C$($?(B. 2 : boundary. 3 : $B$H$J$j$O0[$J$k(B.
3168: % [ ] join $B$r$d$C$F(B $B%Y%/%H%k$N(B clone $B$r:n$k(B.
3169: cone.fan k get , (flipped) getNode 2 get /vv set
3170: /newmerged [ ] vv join def
3171: 0 1 , newmerged length 1 sub {
3172: /i set
3173: newmerged i get , (2).. eq { }
3174: { newmerged i (0).. put } ifelse
3175: } for
3176: % newconeid $B$O(B k $B$G$&$a$k(B.
3177: /newconeid newfacets length newVector { pop k to_univNum } map def
3178:
3179: % merged $B$H(B newmerged $B$r(B cone $B$NNY@\4X78$N$_$G99?7$9$k(B.
3180: % $BF1$8(B init $B$r;}$D$3$H$O$o$+$C$F$$$k$N$G(B facet vector $B$N$_$N(B check $B$G==J,(B.
3181: % merged $B$N(B i $BHVL\(B $B$H(B newmerged $B$N(B j $BHVL\$GHf3S(B.
3182: 0 1 , merged length 1 sub {
3183: /i set
3184: 0 1 , newmerged length 1 sub {
3185: /j set
3186: merged i get , (0).. eq ,
3187: newmerged j get , (0).. eq , and
3188: nextcid i get , k to_univNum eq , and
3189: {
3190: facets i get , newfacets j get , add isZero {
3191: % merged[i], newmerged[j] $B$K(B 1 $B$rF~$l$F>C$9(B.
3192: % $B>e$NH=Dj$O(B nextfid, newnextfid $B$rMQ$$$F$b$h$$$N$G$O(B?
3193: merged i (1).. put
3194: newmerged j (1).. put
3195: } { } ifelse
3196: } { } ifelse
3197: } for
3198: } for
3199:
3200: % Step2. $B7k9g$7$F$+$i(B, $B$^$@D4$Y$F$J$$(B facet $B$rC5$9(B.
3201: cone.debug { (Step 2. Joining *** and new***) message } { } ifelse
3202: /facets facets newfacets join def
3203: /merged merged newmerged join def
3204: /nextcid nextcid newnextcid join def
3205: /nextfid nextfid newnextfid join
3206: /coneid coneid newconeid join def
3207:
3208: cone.debug{ ( Checking facets.) message } { } ifelse
3209: /k null def
3210: 0 1 , merged length 1 sub {
3211: /i set
3212: % i message
3213: merged i get (0).. eq {
3214: % i $BHVL\$r$^$@D4$Y$F$$$J$$(B.
3215: coneid i get , /p set
3216: nextcid i get , /q set
3217: cone.debug { [p q] message } { } ifelse
3218: q (0).. ge {
3219: % cone.fan [p] $B$H(B cone.fan [q] $B$N(B initial $B$rHf3S$9$k(B.
3220: % $BF1$8$J$i(B k $B$r@_Dj(B. exit for. $B0c$($P(B merged[i] = 3 ($B0c$&(B) $B$rBeF~(B.
3221: % differentC $B$O$9$G$K(B $B8=:_$N(B dhcone $B$H0c$&$H(B check $B$5$l$?(B cone $BHV9f(B.
3222: % dhcone.checked $B$O(B dhcone $B$,$9$G$K@8@.$5$l$F$$$k(B cone $BHV9f$N%j%9%H(B.
3223: % $B$3$l$K$O$$$C$F$$$F$b0c$&(B.
3224: q differentC memberQ , q dhcone.checked memberQ , or
3225: { /rr [0 ] def }
3226: { p q isSameCone_h /rr set } ifelse
3227:
3228: rr 0 get 1 eq {
3229: cone.debug { (Found next cone. ) message } { } ifelse
3230: /k q to_univNum def exit
3231: } {
3232: cone.debug { ( It is a different cone. ) message } { } ifelse
3233: differentC [ q ] join /differentC set
3234: merged i (3).. put
3235: } ifelse
3236: } { } ifelse
3237: } { } ifelse
3238: } for
3239:
3240: k tag 0 eq { exit } { } ifelse
3241: } loop
3242:
3243: [(-1)..] cones join shell rest /cones set
3244: % dhcone.checked, dhcone.unchecked $B$r99?7(B.
3245: dhcone.checked cones join /dhcone.checked set
3246: dhcone.unchecked cones setMinus /dhcone.unchecked set
3247:
3248: [(dhcone) [ ]
3249: [
3250: [(cones) [ ] cones] arrayToTree
3251: [(facets) [ ] facets] arrayToTree
3252: [(merged) [ ] merged] arrayToTree
1.5 takayama 3253: [(nextcid) [ ] nextcid] arrayToTree
3254: [(nextfid) [ ] nextfid] arrayToTree
3255: [(coneid) [ ] coneid] arrayToTree
1.4 takayama 3256: ]
3257: ] arrayToTree /arg1 set
3258: ] pop
3259: popVariables
3260: arg1
3261: } def
3262:
3263:
3264: %<
3265: % Usages: dhCones_h
3266: % cone.fan $B$O(B doubly homogenized (local) $B$G@8@.$5$l$?(B Grobner fan.
3267: % cone.fan $B$r(B dehomogenize (H->1) $B$7$F(B init $B$rHf$Y$F(B dhcone.fan $B$r@8@.$9$k(B.
3268: %
3269: % $B%F%9%H(B1. (data/test14.sm1) run (data/test14-data.sm1) run
3270: % dhCones_h
3271: % test22
3272: %>
3273: /dhCones_h {
3274: (cone.grGblist) boundp { } {setGrGblist} ifelse
3275: genDhcone.init
3276: /dhcone.fan [ ] def
3277: {
3278: (-----------------------------------------) message
3279: (#dhcone.unchecked = ) messagen dhcone.unchecked length message
3280: dhcone.unchecked length 0 eq { exit } { } ifelse
3281: dhcone.fan
3282: [ dhcone.unchecked 0 get , genDhcone ] join /dhcone.fan set
3283: (#dhcone.fan = ) messagen dhcone.fan length message
3284: } loop
3285: dhcone.fan
3286: } def
3287:
1.5 takayama 3288: %<
3289: % Usages: dhcone.rtable
3290: % dhcone $B$NHV9f$H(B cone $B$NHV9f$N(B $BCV49I=$r@8@.$7(B dhcone2.fan (merge $B$7$?(B cone $B$N>pJs(B)
3291: % $B$r(B dhcone.fan $B$+$i:n$k(B. dhcone2.gblist $B$b:n$kJd=u4X?t(B.
3292: % dhCones_h $B$7$F$+$i(B dhcone.rable $B$9$k(B.
3293: %>
3294: /dhcone.rtable {
3295: [/i /j /vv /cones /facets /facets2 /merged /nextcid /nextcid2 /ii /ww] pushVariables
3296: [
3297: % $BCV49I=(B dhcone.h2dh $B$r:n$k(B.
3298: /dhcone.h2dh cone.fan length newVector.with-1 def
3299: 0 1 , dhcone.fan length 1 sub {
3300: /i set
3301: dhcone.fan i get , (cones) getNode 2 get /vv set
3302: 0 1 vv length 1 sub {
3303: /j set
3304: dhcone.h2dh , vv j get , i to_univNum , put
3305: } for
3306: } for
3307: % merge $B$7$?(B dhcone $B$r@0M}$7$?$b$N(B, dhcone2.fan $B$r:n$k(B.
3308: /dhcone2.fan dhcone.fan length newVector def
3309: 0 1 , dhcone.fan length 1 sub {
3310: /i set
3311: dhcone.fan i get (facets) getNode 2 get /facets set
3312: dhcone.fan i get (merged) getNode 2 get /merged set
3313: dhcone.fan i get (nextcid) getNode 2 get /nextcid set
3314: dhcone.fan i get (cones) getNode 2 get /cones set
3315: /facets2 [ ] def
3316: /nextcid2 [ ] def
3317: 0 1 , facets length 1 sub {
3318: /j set
3319: merged j get , (3).. eq {
3320: facets2 [ facets j get ] join /facets2 set
3321: % $B$H$J$j$N(B cone $B$,$"$k$H$-(B $BJQ49I=$K$7$?$,$$(B, cone $BHV9f$rJQ49(B
3322: nextcid2 [ dhcone.h2dh , nextcid j get , get ] join /nextcid2 set
3323: } { } ifelse
3324: merged j get , (2).. eq {
3325: facets2 [ facets j get ] join /facets2 set
3326: % $B6-3&$N$H$-(B -2 $B$rF~$l$k(B.
3327: nextcid2 [ (-2).. ] join /nextcid2 set
3328: } { } ifelse
3329: } for
3330:
3331: dhcone2.fan i ,
3332: [(dhcone) [ ]
3333: [
3334: [(facets) [ ] facets2] arrayToTree
3335: [(nextcid) [ ] nextcid2] arrayToTree
3336: [(cones) [ ] cones] arrayToTree
3337: ]
3338: ] arrayToTree , put
3339:
3340: } for
3341:
3342: % $B:G8e$K(B dhcone2.gblist $B$r:n$k(B.
3343: /dhcone2.gblist , dhcone2.fan length newVector , def
3344: 0 1 , dhcone2.fan length 1 sub {
3345: /i set
3346: dhcone2.fan i get (cones) getNode 2 get /cones set
3347: cone.grGblist , cones 0 get , get , /ii set % GB of initial (H->1).
3348: cone.gblist i get , (weight) getNode , [ 2 0 2 ] get /ww set
3349:
3350: dhcone2.gblist i,
3351: [(gbasis) [ ]
3352: [
3353: [(initial) [ ] ii] arrayToTree
3354: [(weight) [ ] ww] arrayToTree
3355: ]
3356: ] arrayToTree , put
3357:
3358: } for
3359: (dhcone2.fan, dhcone2.gblist, dhcone.h2dh are set.) message
3360:
3361: ] pop
3362: popVariables
3363: } def
3364:
3365: %<
3366: % $BI=$N8+J}$N2r@b$r0u:~$9$k4X?t(B.
3367: % Usages: dhcone.explain
3368: %>
3369: /dhcone.explain {
3370: [
3371: ( ) nl
3372: (Data format in << dhcone2.fan >>, which is a dehomogenized Grobner fan.) nl nl
3373: (<< cone.vlist >> is the list of the variables.) nl
3374: @@@.Hsymbol ( is the homogenization variable to be dehomogenized.) nl nl
3375: (<< cone.input >> is generators of a given ideal.) nl nl
3376: (<< cone.d >> is the dimension of parametrization space of the weights P_w) nl
3377: ( P_w is a cone in R^m where the number m is stored in << cone.m >>) nl
3378: ( P_w --- W ---> R^n [weight space]. ) nl
3379: ( W is stored in << cone.W >> ) nl
3380: ( << u cone.W mul >> gives the weight vector standing for u) nl nl
3381: (All cones in the data lie in the weight parametrization space P_w.) nl
3382: ( "facets" are the inner normal vector of the cone. ) nl
3383: ( "nextcid" is a list of the cone id's of the adjacent cones.) nl
3384: ( -2 in "nextcid" means that this facet lies on the border of the weight space.) nl
3385: ( "cones" is a list of the cone id's of the NON-dehomonized Grobner fan) nl
3386: ( stored in << cone.fan >>) nl
3387: ] cat
3388: } def
3389:
3390: %<
3391: % dhcone.printGrobnerFan
3392: % dhcone $B$N0u:~4X?t(B
3393: %>
3394: /dhcone.printGrobnerFan {
3395: [/i] pushVariables
3396: [
3397: (========== Grobner Fan (for dehomogenized cones) ============) message
3398: [
3399: (cone.comment)
3400: (cone.vlist) (cone.vv)
3401: (cone.input)
3402: (cone.type) (cone.local) (cone.h0)
3403: (cone.n) (cone.m) (cone.d)
3404: (cone.W) (cone.Wpos) (cone.Wt)
3405: (cone.L) (cone.Lp) (cone.Lpt)
3406: (cone.weightBorder)
3407: (cone.incidence)
3408: ] { printGrobnerFan.1 } map
3409: ( ) message
3410: (The number of cones = ) messagen dhcone.fan length message
3411: ( ) message
3412: 0 1 dhcone2.fan length 1 sub {
3413: /ii set
3414: ii messagen ( : ) messagen
3415: dhcone2.fan ii get printTree
3416: } for
3417: 1 {
3418: 0 1 dhcone2.gblist length 1 sub {
3419: /ii set
3420: ii messagen ( : ) messagen
3421: dhcone2.gblist ii get printTree
3422: } for
3423: } { } ifelse
3424:
3425:
3426: (=========================================) message
3427: %(cone.withGblist = ) messagen cone.withGblist message
3428: dhcone.explain message
3429: ( ) message
3430: ] pop
3431: popVariables
3432: } def
3433:
3434: %
3435: % $B;n$7J}(B test14, 22, 25
3436: %
3437: % (data/test14.sm1) run (data/test14-data.sm1) run
3438: % printGrobnerFan ; % H $BIU$-$G0u:~(B.
3439: % dhCones_h ; % dehomogenize Cones.
3440: % dhcone.rtable ; % dhcone2.fan $BEy$r@8@.(B.
3441: % dhcone.printGrobnerFan ; % $B0u:~(B.
3442: % $B0u:~$7$?$b$N$O(B test*-print.txt $B$X3JG<$7$F$"$k(B.
3443: %
3444:
3445: % Todo: save functions.
1.9 takayama 3446:
3447: %<
3448: % Collart, Kalkbrener, Mall $B$N%"%k%4%j%:%`$K$h$k(B gb $B$N(B flip.
3449: % See also Sturmfels' book, p.22, 23.
3450: % Usages: [reducedGb, vlist, oldWeight, facetWeight, newWeight] ckmFlip rGb
3451: % If it fails, then it returns null, else it returns the reducedGb for the
3452: % newWeight.
3453: % gb $B$N(B check $B$r$d$k$N$G(B, $B$=$l$K<:GT$7$?$i(B null $B$rLa$9(B.
3454: % weight $B$O$9$Y$F(B vw $B7A<0$G(B. vw $B7A<0(B = variable weight $B$N7+$jJV$7$N7A<0(B
3455: % reducedGb $B$OJ8;zNs$N%j%9%H$G$O$J$/B?9`<0$N7A<0$N$3$H(B.
3456: %>
3457: /ckmFlip {
3458: /arg1 set
3459: [/arg_ckmFlip /gOld /vlist /oldWeight /facetWeight /newWeight
3460: /gNew
3461: /ww /ww1 /ww2 % $BK\$NCf$N(B w1, w, w2 ($B8E$$(B, facet, $B?7$7$$(B)
3462: /ch1 /ch2 % $BK\$NCf$N(B {\cal H}_1, {\cal H}_2
3463: /grData /rTable
3464: /rTable2 % rTable $B$NH?BP$NJQ49(B.
3465: /facetWeight_gr /vlist_gr % graded ring $BMQ(B.
3466: /oldWeight_gr
3467: /ccf % reduction $B$7$?78?t(B.
3468: /rwork /ccf2 /gNew
3469: ] pushVariables
3470: [
3471: arg1 /arg_ckmFlip set
3472: arg_ckmFlip 0 get /gOld set
3473: arg_ckmFlip 1 get /vlist set
3474: arg_ckmFlip 2 get /oldWeight set
3475: arg_ckmFlip 3 get /facetWeight set
3476: arg_ckmFlip 4 get /newWeight set
3477:
3478: % facet weight vector ww $B$K$D$$$F$N(B initial $B$r<h$j=P$9(B. ch1 $B$X$$$l$k(B.
3479: gOld getRing ring_def
3480: facetWeight weightv /ww set
3481: gOld { ww init } map /ch1 set % facetWeight $B$K$h$k(B initial $B$N<h$j=P$7(B.
3482:
3483:
3484: % $BNc(B: [(x,y) [(x) -1 (Dx) 1 (y) -1 (Dy) 2]] getGrRing
3485: % [$x,y,y',$ , [ $x$ , $y$ ] , [ [ $Dy$ , $y'$ ] ] ]
3486: % $BJQ?t%j%9%H(B $BCV49I=(B
3487: % ch1 $B$r(B gr_ww $B$N85$KJQ49(B.
3488: [vlist facetWeight] getGrRing /grData set
3489: [grData 0 get ring_of_differential_operators 0] define_ring /rwork set
3490: grData 2 get { { . } map } map /rTable set
3491: rTable { reverse } map /rTable2 set
3492: grData 0 get /vlist_gr set
3493: ch1 { toString . rTable replace toString } map /ch1 set
3494:
3495: oldWeight { dup isString { . rTable replace toString }
3496: { } ifelse } map /oldWeight_gr set
3497:
3498: % facetWeight $B$b(B $B?7$7$$4D(B gr_ww $B$N(B weight $B$KJQ49(B.
3499: % $BNc(B. [(x) -1 (Dx) 1 (y) -1 (Dy) 2] ==> [(x) -1 (Dx) 1 (y) -1 (y') 2]
3500: facetWeight { dup isString { . rTable replace toString }
3501: { } ifelse } map /facetWeight_gr set
3502: % Dx x = x Dx + h H or Dx x = x Dx + h^2 $B$G7W;;(B.
3503: % $B$I$A$i$r$H$k$+$O(B cone.gb_gr $B$G6hJL$9$k$7$+$J$7(B
3504: %% [ch1 vlist_gr oldWeight_gr] /ttt set
3505: %% ttt cone.gb_gr /ch1 set %$B:FEY$N7W;;$OITMW(B.
3506: [[(1)] vlist_gr oldWeight_gr] cone.gb_gr getRing ring_def % Set Ring.
3507: ch1 {toString .} map /ch1 set
3508: %% $B$3$3$^$G$G$H$j$"$($:%F%9%H$r$7$h$&(B.
3509: %% ch1 /arg1 set
3510: % newWeight $B$b(B $B?7$7$$4D(B gr_ww $B$N(B weight $B$KJQ49(B.
3511: % $BNc(B. [(x) -1 (Dx) 1 (y) -1 (Dy) 2] ==> [(x) -1 (Dx) 1 (y) -1 (y') 2]
3512: newWeight { dup isString { . rTable replace toString }
3513: { } ifelse } map /newWeight_gr set
3514: [ch1 { toString } map vlist_gr newWeight_gr] cone.gb_gr /ch2 set
3515:
3516: % Dx x = x Dx + h H or Dx x = x Dx + h^2 $B$G7W;;(B.
3517: % $B$I$A$i$r$H$k$+$O(B cone.reduction_gr $B$G6hJL$9$k$7$+$J$7(B
3518: ch1 getRing ring_def ;
3519: ch2 {toString .} map {ch1 cone.reduction} map /ccf set
3520: %ccf pmat
3521: % $B$H$j$"$($:%F%9%H(B.
3522: % [ch1 ch2] /arg1 set
3523: %% ccf[i][0] $B$O(B 0 $B$G$J$$$HL7=b(B. check $B$^$@$7$F$J$$(B.
3524:
3525: %% ccf[i][2] (syzygy) $B$r(B gr $B$+$i(B $B$b$H$N(B ring $B$XLa$7(B,
3526: %% $B?7$7$$(B reduced gbasis $B$r(B ccf[i][2] * gOld $B$G:n$k(B.
3527: rwork ring_def
3528: ccf { 2 get {toString . rTable2 replace toString} map } map /ccf2 set
3529: %% ccf2 $B$O(B gr $B$G$J$$(B ring $B$N85(B.
3530: gOld getRing ring_def
1.10 ! takayama 3531: cone.DhH { cone.begin_DhH } { } ifelse % Hh $B$+(B h^2 $B$+(B.
1.9 takayama 3532: ccf2 { {.} map gOld mul } map /gNew set
3533: gNew { toString } map /gNew set
1.10 ! takayama 3534: cone.DhH { cone.end_DhH } { } ifelse % Hh $B$+(B h^2 $B$+(B.
1.9 takayama 3535: % gNew /arg1 set
3536: %gNew $B$,(B newWeight $B$G$N(B GB $B$+(B check. Yes $B$J$i(B reduced basis $B$X(B.
3537: %No $B$J$i(B null $B$rLa$9(B.
1.10 ! takayama 3538: %%Ref: note @s/2005/06/30-note-gfan.pdf
1.9 takayama 3539: gNew [(gbCheck) 1] setAttributeList newWeight
3540: cone.gb (gb) getAttribute
3541: 1 eq {
3542: gNew [(reduceOnly) 1] setAttributeList newWeight cone.gb /arg1 set
3543: }{ /arg1 null def } ifelse
3544: ] pop
3545: popVariables
3546: arg1
3547: } def
3548:
3549: %<
3550: % Usages: f gbasis cone.reduction_DhH
1.10 ! takayama 3551: % dx x = x dx + h H $B$G$N(B reduction.
1.9 takayama 3552: %>
3553: /cone.reduction_DhH {
3554: /arg2 set /arg1 set
3555: [/ff /ggbasis /eenv /ans] pushVariables
3556: [
3557: /ff arg1 def /ggbasis arg2 def
1.10 ! takayama 3558: cone.begin_DhH
! 3559: ff ggbasis reduction /ans set
! 3560: cone.end_DhH
! 3561: /arg1 ans def
! 3562: ] pop
! 3563: popVariables
! 3564: arg1
! 3565: } def
! 3566:
! 3567: %<
! 3568: % Usages: f gbasis cone.reduction_Dh
! 3569: % dx x = x dx + h^2 $B$G$N(B reduction.
! 3570: %>
! 3571: /cone.reduction_Dh {
! 3572: /arg2 set /arg1 set
! 3573: [/ff /ggbasis /eenv /ans] pushVariables
! 3574: [
! 3575: /ff arg1 def /ggbasis arg2 def
1.9 takayama 3576: ff ggbasis reduction /ans set
3577: /arg1 ans def
3578: ] pop
3579: popVariables
3580: arg1
3581: } def
3582:
1.10 ! takayama 3583: %<
! 3584: % Usages: cone.begin_DhH dx x = x dx + h H $B$r3+;O(B.
! 3585: %>
1.9 takayama 3586: /cone.begin_DhH {
3587: [(Homogenize) (AutoReduce) (KanGBmessage)] pushEnv /cone.eenv set
3588: [(Homogenize) 3] system_variable
3589: } def
3590:
1.10 ! takayama 3591: %<
! 3592: % Usages: cone.begin_DhH dx x = x dx + h H $B$r=*N;(B.
! 3593: %>
1.9 takayama 3594: /cone.end_DhH {
3595: cone.eenv popEnv
3596: } def
3597:
1.10 ! takayama 3598: %<
! 3599: % Usages: ff vv ww cone.gb_gr_DhH dx x = x dx + h H $B$G7W;;(B.
! 3600: % dh.gb $B$O(B dhecart.sm1 $B$GDj5A$5$l$F$*$j(B, dx x = x dx + h H $B$G$N7W;;(B.
! 3601: % gr $B$r$H$C$F$b(B, -w,w $B$N>l9g$O(B $BHyJ,:nMQAG4D$N$^$^$G$"$j(B, $B$3$l$,I,MW(B.
! 3602: % bug? cone.gb $B$G==J,(B?
! 3603: %>
! 3604: /cone.gb_gr_DhH {
! 3605: /arg1 set
! 3606: [/ff /ww /vv] pushVariables
! 3607: [
! 3608: /ff arg1 0 get def
! 3609: /vv arg1 1 get def
! 3610: /ww arg1 2 get def
! 3611: /dh.gb.verbose 1 def
! 3612: /dh.autoHomogenize 0 def
! 3613: [(AutoReduce) 1] system_variable
! 3614: [ff { toString } map vv
! 3615: [ww vv generateD1_1]] dh.gb 0 get /arg1 set
! 3616: ] pop
! 3617: popVariables
! 3618: arg1
! 3619: } def
! 3620: %<
! 3621: % Usages: ff vv ww cone.gb_gr_Dh dx x = x dx + h^2 $B$G7W;;(B.
! 3622: % gb $B$O(B dhecart.sm1 $B$GDj5A$5$l$F$*$j(B, dx x = x dx + h^2 $B$G$N7W;;(B.
! 3623: % gr $B$r$H$C$F$b(B, -w,w $B$N>l9g$O(B $BHyJ,:nMQAG4D$N$^$^$G$"$j(B, $B$3$l$,I,MW(B.
! 3624: % bug? cone.gb $B$G==J,(B?
! 3625: %>
! 3626: /cone.gb_gr_Dh {
! 3627: /arg1 set
! 3628: [/ff /ww /vv] pushVariables
! 3629: [
! 3630: /ff arg1 0 get def
! 3631: /vv arg1 1 get def
! 3632: /ww arg1 2 get def
! 3633: /gb.verbose 1 def
! 3634: /gb.autoHomogenize 0 def
! 3635: [(AutoReduce) 1] system_variable
! 3636: [ff { toString } map vv
! 3637: [ww vv generateD1_1]] gb 0 get /arg1 set
! 3638: /gb.autoHomogenize 1 def
! 3639: ] pop
! 3640: popVariables
! 3641: arg1
! 3642: } def
! 3643:
! 3644:
! 3645: % $B$3$l$i$O(B cone.ckmFlip 1 $B$N;~$7$+;H$o$:(B.
! 3646: /cone.reduction {
! 3647: cone.DhH {
! 3648: cone.reduction_DhH
! 3649: }{
! 3650: cone.reduction_Dh
! 3651: } ifelse
! 3652: } def
! 3653: /cone.gb_gr {
! 3654: cone.DhH {
! 3655: cone.gb_gr_DhH
! 3656: }{
! 3657: cone.gb_gr_Dh
! 3658: } ifelse
! 3659: } def
! 3660:
! 3661:
1.9 takayama 3662: /test1.ckmFlip {
3663: % cf. cone.sample2
3664: cone.load.cohom
3665: /cone.comment [
3666: (BS for y and y-(x-1)^2, t1, t2 space, in doubly homogenized Weyl algebra.) nl
3667: (The Grobner cones are dehomogenized to get local Grobner fan.) nl
3668: ] cat def
3669: /cone.vlist [(t1) (t2) (x) (y) (Dt1) (Dt2) (Dx) (Dy) (h) (H)] def
3670: /cone.vv (t1,t2,x,y) def
3671: /cone.type 1 def
3672: /cone.parametrizeWeightSpace {
3673: 4 2 parametrizeSmallFan
3674: } def
1.10 ! takayama 3675:
! 3676: /cone.DhH 1 def
! 3677: /cone.ckmFlip 1 def
! 3678:
1.9 takayama 3679: /cone.local 1 def
3680: /cone.w_start null def
3681: /cone.h0 1 def
3682: /cone.input
3683: [
3684: (t1-y) (t2 - (y-(x-1)^2))
3685: ((-2 x + 2)*Dt2+Dx)
3686: (Dt1+Dt2+Dy)
3687: ]
3688: def
3689: % homogenize
3690: [cone.vv ring_of_differential_operators
3691: [[(t1) -1 (t2) -1 (Dt1) 1 (Dt2) 1]] ecart.weight_vector
3692: 0] define_ring
3693: dh.begin
3694: cone.input { . homogenize toString } map /cone.input set
3695: dh.end
3696:
3697:
3698: % $B%F%9%H$r3+;O$9$k(B.
3699: % getStartingCone /cone.ncone set
3700: % cone.ncone updateFan
3701: % cone.gblist 0 get message
3702: % cone.ncone /cone.ccone set
3703: % getNextFlip /cone.nextflip set
3704: % cone.nextflip message
3705:
3706: /wOld [(t1) , -29 , (t2) , -38 , (Dt1) , 29 , (Dt2) , 38 ] def
3707: /wFacet [(t1) , -1 , (t2) , -1 , (Dt1) , 1 , (Dt2) , 1 ] def
3708: /wNew [(t1) , -39 , (t2) , -38 , (Dt1) , 39 , (Dt2) , 38 ] def
3709: cone.input wOld cone.gb /ff set
3710: [ff (t1,t2,x,y) wOld wFacet wNew] ckmFlip /ff2 set
3711: (See ff and ff2) message
3712:
3713: } def
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>