=================================================================== RCS file: /home/cvs/OpenXM/src/kan96xx/Doc/gfan.sm1,v retrieving revision 1.10 retrieving revision 1.20 diff -u -p -r1.10 -r1.20 --- OpenXM/src/kan96xx/Doc/gfan.sm1 2005/07/07 01:31:21 1.10 +++ OpenXM/src/kan96xx/Doc/gfan.sm1 2018/05/02 02:28:13 1.20 @@ -1,9 +1,43 @@ -% $OpenXM: OpenXM/src/kan96xx/Doc/gfan.sm1,v 1.9 2005/06/30 08:39:39 takayama Exp $ +% $OpenXM: OpenXM/src/kan96xx/Doc/gfan.sm1,v 1.19 2013/10/11 01:08:35 takayama Exp $ % cp cone.sm1 $OpenXM_HOME/src/kan96xx/Doc/gfan.sm1 -% $Id: gfan.sm1,v 1.10 2005/07/07 01:31:21 takayama Exp $ +% $Id: gfan.sm1,v 1.20 2018/05/02 02:28:13 takayama Exp $ % iso-2022-jp %%Ref: @s/2004/08/21-note.pdf +%% gfan.sm1 works only for polymake 2.0 Use webservice of 2.0. +[(gfan) +[ + (gfan.sm1 is a package to compute global and local Grobner fans.) + (See R.Bahloul and N.Takayama, arxiv, math.AG/0412044 and references as to algorithms.) + (At the beginning of the source code gfan.sm1, there are sample inputs cone.sample and cone.sample2.) + ( ) + (gfan.sm1 works only with polymake 2.0. We provide a web service of computing ) + (with polymake 2.0. /@@@polymake.web 1 def is set by default in gfan.sm1.) + (See changelog-ja.tex as to details on the difference between 2.0 and later versions.) + ( ) + (*cone.sample ; is an example. See the source code. The state polytope is the hexagon.) + ( ) + (*cone.Wt cone.Lpt {vertices in the output} are weights on the rays of the Grobner cone.) + (*cone.L gives a basis of the linearity space.) + (*cone.Lp gives a basis of the pointed cone. cone.Lpt is the transpose of cone.Lp.) + $*When v is a row vector in an ouput cone, (v cone.Lp cone.W) gives $ + ( the corresponding weight vector in the full variable space in D) + (*cone.incidence is a list of [[cone num1,facet num1], [cone num2,facet num2]]) + ( which means that cone num1 and cone num2 are adjacent and shares ) + ( the facet num1 and the facet num2) + (*/cone.withGblist 1 def saves the Grobner basis standing for each cone.) + ( ) + (*Cone descriptions: cone.fan) + (**A facet is given by its normal vector n of a cone. It gives facet num of the cone.) + (**A cone is defined by facet normal vectors n1, n2, ... as n1.x>=0 and n2.x >=0 and ...) + (**facetsv is a list of facets expressed by generators.) + (**nextcid is a list of the adjacent cone numbers.) + (**nextfid is a list of the shared facet numbers.) + (**vertices is the generators of a cone.) + (**inequalities are not necessarily unique.) +] +] putUsages + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% Two examples are given below to get a global Grobner fan and %% a local Grobner fan ; cone.sample and cone.sample2 @@ -17,6 +51,7 @@ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% /cone.sample { cone.load.cohom + /cone.ckmFlip 1 def % write a comment about the problem. "nl" means new line. /cone.comment [ (Toric ideal for 1-simplex x 2-simplex, in k[x]) nl @@ -90,6 +125,9 @@ def cone.comment message (cone.input = ) message cone.input message +%%% Step 0. If you want to output Grobner basis standing for each cone, then uncomment +% /cone.withGblist 1 def + %%%% Step 1. Enumerating the Grobner Cones in a global ring. %%%% The result is stored in cone.fan getGrobnerFan @@ -98,7 +136,7 @@ getGrobnerFan printGrobnerFan %%%% If you want to save the data to the file sm1out.txt, then uncomment. -% /cone.withGblist 1 def saveGrobnerFan /ff set ff output +%saveGrobnerFan /ff set ff output %%%% Step 2. Dehomogenize the Grobner Cones %%%% by the equivalence relation in a local ring (uncomment). @@ -129,6 +167,7 @@ printGrobnerFan %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% /cone.sample2 { cone.load.cohom + /cone.ckmFlip 1 def % write a comment about the problem. "nl" means new line. /cone.comment [ (BS for y and y-(x-1)^2, t1, t2 space, in doubly homogenized Weyl algebra.) nl @@ -248,14 +287,16 @@ dhcone.printGrobnerFan % If you use the cgi/polymake on the net, then uncomment out. %/doPolymake {doPolymake.OoHG} def (Using doPolymake.OoHG ) message %/polymake.start {polymake.start.OoHG} def (Using polymake.start.OoHG ) message +/@@@polymake.web 1 def %% Choose it automatically. -[(which) (polymake)] oxshell tag 0 eq { - (Polymake is not installed in this system.) message - /doPolymake {doPolymake.OoHG} def - (Using doPolymake.OoHG ) message - /polymake.start {polymake.start.OoHG} def - (Using polymake.start.OoHG ) message -} { (Local polymake will be used.) message } ifelse +[(which) (polymake)] oxshell tag 0 eq +@@@polymake.web 1 eq +or +{ + (Polymake is not installed in this system or @@@polymake.web is set.) message + usePolymake.OoHG.curl + (Using doPolymake.OoHG.curl ) message +} { usePolymake.local (Local polymake will be used.) message } ifelse /cone.debug 1 def @@ -413,6 +454,13 @@ cone.comment message %> /cone.DhH 0 def +%< +% Global +% gbCheck をするか? しないと結果はあやふや. しかしメモリ exhaust は防げる. +% 使うときは /cone.epsilon, /cone.epsilon.limit を十分小さくしておく. +%> +/cone.do_gbCheck 1 def + % Default の cone.gb の定義. 各プログラムで再度定義してもよい. /cone.gb { cone.DhH { @@ -1011,17 +1059,31 @@ def %> /getConeInfo { /arg1 set - [/ww /g /ceq /ceq2 /cdim /mmc /mmL /rr /ineq /ppt] pushVariables + [/ww /g /ceq /ceq2 /cdim /mmc /mmL /rr /ineq /ppt /rr0 /mm0 /mm1] pushVariables [ /ceq arg1 def ceq pruneZeroVector /ceq set + + ceq length 0 eq { + (Monomial ideal is not accepted as an input.) cone_ir_input + } { } ifelse + + /mm1 + ( Use [(keep_tmp_files) 1] oxshell to check the input to polymake2tfb. See /tmp or $TMP ) + def + ceq genPo2 /ceq2 set % ceq2 は polymake.data(polymake.INEQUALITIES(...)) 形式 % polymake で ceq2 の次元の計算. /getConeInfo.ceq ceq def /getConeInfo.ceq2 ceq2 def cone.debug { (Calling polymake DIM.) message } { } ifelse - [(DIM) ceq2] doPolymake 1 get /rr set + [(DIM) ceq2] doPolymake /rr0 set + % rr0 2 get message + rr0 2 get 1 get 0 get /mm0 set + mm0 length 0 eq { } + { [mm0 mm1] cat error } ifelse + rr0 1 get /rr set cone.debug {(Done.) message } { } ifelse % test5 には次のコメントとりさる. 上の行をコメントアウト. % test5.data tfbToTree /rr set @@ -1036,13 +1098,24 @@ def % FACETS を持っていないなら再度計算する. % POINTED, NOT__POINTED も得られる cone.debug { (Calling polymake FACETS.) message } { } ifelse - [(FACETS) ceq2] doPolymake 1 get /rr set + [(FACETS) ceq2] doPolymake /rr0 set + + % rr0 2 get message + rr0 2 get 1 get 0 get /mm0 set + mm0 length 0 eq { } + { [mm0 mm1] cat error } ifelse + + rr0 1 get /rr set cone.debug { (Done.) message } { } ifelse } { } ifelse rr (VERTICES) getNode tag 0 eq { (internal error: VERTICES is not found.) error - } { } ifelse + } { + rr (VERTICES) getNode + (UNDEF) getNode tag 0 eq { } + { (internal error: VERTICES is UNDEF. See rr. Set /@@@polymake.web 1 def) error } ifelse + } ifelse /cone.getConeInfo.rr1 rr def @@ -1612,7 +1685,7 @@ def /vlist arg1 def /wlist arg2 def wlist length vlist length eq { - } { (cone_wtowv: length of the argument must be the same.) error} ifelse + } { (cone_wtowv: length of the argument must be the same. Please check the values of cone.vlist cone.vv cone.type parametrizeWeightSpace) error} ifelse wlist to_int32 /wlist set [ @@ -1869,6 +1942,11 @@ def % note: 2004.9.2 cone (facetsv) getNode 2 get facet_i get /v set cone (facets) getNode 2 get facet_i get /f set + + v length 0 eq { + (The codimension of the linarity space of the Grobner cone seems to be 1 or 0.) cone_ir_input + } { } ifelse + /vp v 0 get def 1 1 v length 1 sub { /ii set @@ -2046,15 +2124,20 @@ def %> /cone.gb_Dh { /arg2 set /arg1 set - [/ff /ww /gg] pushVariables + [/ff /ww /gg /gbopt] pushVariables [ /ff arg1 def /ww arg2 def [(AutoReduce) 1] system_variable [cone.vv ring_of_differential_operators [ww] weight_vector 0] define_ring - [ff {toString .} map] ff getAttributeList setAttributeList - groebner 0 get /gg set + %(---) messagen ff getAttributeList message + ff getAttributeList tag 0 eq {/gbopt [ ] def } + { + /gbopt ff getAttributeList def + } ifelse + [ff {toString .} map gbopt] + groebner 0 get /gg set %% groenber は attribute を受け付けない. /cone.gb_Dh.g gg def /arg1 gg def ] pop @@ -2384,7 +2467,7 @@ def % Usages: result_getNextFlip getNextCone ncone % flip して新しい ncone を得る. %> -/getNextCone { +/getNextCone.orig { /arg1 set [/ncone /ccone /kk /w /next_weight_w_wv] pushVariables [ @@ -2718,6 +2801,7 @@ def /printGrobnerFan { [/i] pushVariables [ + $(gfan) usage to find explanations on variables.$ message (========== Grobner Fan ====================) message [ (cone.comment) @@ -3453,6 +3537,7 @@ def % gb の check をやるので, それに失敗したら null を戻す. % weight はすべて vw 形式で. vw 形式 = variable weight の繰り返しの形式 % reducedGb は文字列のリストではなく多項式の形式のこと. +% 理由は reducedGb より ring の構造を読むため. %> /ckmFlip { /arg1 set @@ -3499,6 +3584,12 @@ def % 例. [(x) -1 (Dx) 1 (y) -1 (Dy) 2] ==> [(x) -1 (Dx) 1 (y) -1 (y') 2] facetWeight { dup isString { . rTable replace toString } { } ifelse } map /facetWeight_gr set + +% newWeight も 新しい環 gr_ww の weight に変換. +% 例. [(x) -1 (Dx) 1 (y) -1 (Dy) 2] ==> [(x) -1 (Dx) 1 (y) -1 (y') 2] + newWeight { dup isString { . rTable replace toString } + { } ifelse } map /newWeight_gr set + % Dx x = x Dx + h H or Dx x = x Dx + h^2 で計算. % どちらをとるかは cone.gb_gr で区別するしかなし %% [ch1 vlist_gr oldWeight_gr] /ttt set @@ -3507,10 +3598,6 @@ def ch1 {toString .} map /ch1 set %% ここまででとりあえずテストをしよう. %% ch1 /arg1 set -% newWeight も 新しい環 gr_ww の weight に変換. -% 例. [(x) -1 (Dx) 1 (y) -1 (Dy) 2] ==> [(x) -1 (Dx) 1 (y) -1 (y') 2] - newWeight { dup isString { . rTable replace toString } - { } ifelse } map /newWeight_gr set [ch1 { toString } map vlist_gr newWeight_gr] cone.gb_gr /ch2 set % Dx x = x Dx + h H or Dx x = x Dx + h^2 で計算. @@ -3536,8 +3623,15 @@ def %gNew が newWeight での GB か check. Yes なら reduced basis へ. %No なら null を戻す. %%Ref: note @s/2005/06/30-note-gfan.pdf - gNew [(gbCheck) 1] setAttributeList newWeight - cone.gb (gb) getAttribute + cone.do_gbCheck not { + (Warning! gbCheck is skipped.) message + } { + (Doing gbCheck.) message + } ifelse + cone.do_gbCheck { + gNew [(gbCheck) 1] setAttributeList newWeight + cone.gb (gb) getAttribute + } { 1 } ifelse 1 eq { gNew [(reduceOnly) 1] setAttributeList newWeight cone.gb /arg1 set }{ /arg1 null def } ifelse @@ -3625,17 +3719,22 @@ def %> /cone.gb_gr_Dh { /arg1 set - [/ff /ww /vv] pushVariables + [/ff /ww /vv /gg /envtmp] pushVariables [ /ff arg1 0 get def /vv arg1 1 get def /ww arg1 2 get def - /gb.verbose 1 def - /gb.autoHomogenize 0 def + + [(AutoReduce) (KanGBmessage)] pushEnv /envtmp set [(AutoReduce) 1] system_variable - [ff { toString } map vv - [ww vv generateD1_1]] gb 0 get /arg1 set - /gb.autoHomogenize 1 def + [(KanGBmessage) 1] system_variable + [vv ring_of_differential_operators + [ww] weight_vector 0] define_ring + [ff {toString .} map] ff getAttributeList setAttributeList + groebner 0 get /gg set + envtmp popEnv + + /arg1 gg def ] pop popVariables arg1 @@ -3710,4 +3809,125 @@ def [ff (t1,t2,x,y) wOld wFacet wNew] ckmFlip /ff2 set (See ff and ff2) message +} def + +%< +% Usages: cone i getaVectorOnFacet +% cone の i 番目の facet の上の vector を求める. +% cf. liftWeight +%> +/getaVectorOnFacet { + /arg2 set /arg1 set + [/cone /facet_i /ep /vp /v /v /ii] pushVariables + [ + /cone arg1 def /facet_i arg2 def + facet_i to_int32 /facet_i set + + cone (facetsv) getNode 2 get facet_i get /v set + /vp v 0 get def + 1 1 v length 1 sub { + /ii set + vp v ii get add /vp set + } for + vp nnormalize_vec /vp set + /arg1 vp def + ] pop + popVariables + arg1 +} def + +/getNextCone { + getNextCone_ckm +} def + +%< +% Usages: result_getNextFlip getNextCone_ckm ncone +% flip して新しい ncone を得る. Collar-Kalkbrener-Moll のアルゴリズムを使う +% if (cone.ckmFlip == 0) 普通の計算 else CKM. +%> +/getNextCone_ckm { + /arg1 set + [/ncone /ccone /kk /w /next_weight_w_wv /cid /ttt] pushVariables + [ + /ccone arg1 def + /ncone null def + /kk ccone 1 get def % kk は cid 番目の cone の kk 番目の facet を表す. + /cid ccone 2 get def % cid は cone の 番号. + ccone 0 get /ccone set + { + ccone tag 0 eq { exit } { } ifelse + +% ccone の kk 番目の facet について flip する. + ccone kk cone.epsilon flipWeight /w set + (Trying new weight is ) messagen w message + w liftWeight /next_weight_w_wv set + (Trying new weight [w,wv] is ) messagen next_weight_w_wv message + + cone.ckmFlip { + [ + cone.gblist cid get (grobnerBasis) getNode 2 get % reduce gb + cone.vv + cone.gblist cid get (weight) getNode [2 0 2] get % weight + ccone kk getaVectorOnFacet liftWeight 1 get % weight on facet + next_weight_w_wv 1 get % new weight + ] /ttt set + ttt message + ttt ckmFlip /cone.cgb set + }{ + cone.input next_weight_w_wv 1 get cone.gb /cone.cgb set + } ifelse + + cone.cgb tag 0 eq not { + [w] next_weight_w_wv join /cone.cgb_weight set + next_weight_w_wv 1 get cone.cgb coneEq /cone.g_ineq set + cone.g_ineq cone.w_ineq join cone.Wt mul cone.Lpt mul + pruneZeroVector /cone.gw_ineq_projectedWtLpt set + + (cone.gw_ineq_projectedWtLpt is obtained.) message + + cone.gw_ineq_projectedWtLpt getConeInfo /cone.nextConeInfo set +% 次元を調べる. だめなら retry + cone.nextConeInfo 0 get 0 get to_int32 cone.d eq { + cone.nextConeInfo 1 get newCone /ncone set + ccone ncone getCommonFacet 0 get { + (Flip succeeded.) message + exit + } { } ifelse + } { } ifelse +% common face がなければ やはり epsilon を小さく. + cone.nextConeInfo 0 get 0 get to_int32 cone.d eq { + (ccone and ncone do not have a common facet.) message + } { + (ncone is not maximal dimensional. ) message + } ifelse + }{ } ifelse + + (Decreasing epsilon to ) messagen + cone.epsilon (1).. (2).. div mul /cone.epsilon set + cone.epsilon cone.epsilon.limit sub numerator (0).. lt { + (Too small cone.epsilon ) error + } { } ifelse + cone.epsilon message + } loop + /arg1 ncone def + ] pop + popVariables + arg1 +} def + +%%change +/cone_ir_input { + /arg1 set + [/msg ] pushVariables + [ + /msg arg1 def + (---------------) message + msg message + ( ) message + (Please also refer to the value of the variables cone.getConeInfo.rr0) message + ( cone.getConeInfo.rr1 cone.Lp cone.cinit) message + $ cone.cinit (FACETS) getNode :: $ message + (We are sorry that we cannot accept this input.) error + ] pop + popVariables } def