=================================================================== RCS file: /home/cvs/OpenXM/src/kan96xx/Doc/gfan.sm1,v retrieving revision 1.7 retrieving revision 1.20 diff -u -p -r1.7 -r1.20 --- OpenXM/src/kan96xx/Doc/gfan.sm1 2004/09/30 07:45:04 1.7 +++ OpenXM/src/kan96xx/Doc/gfan.sm1 2018/05/02 02:28:13 1.20 @@ -1,8 +1,43 @@ -% $OpenXM: OpenXM/src/kan96xx/Doc/gfan.sm1,v 1.6 2004/09/30 07:39:42 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.7 2004/09/30 07:45:04 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 @@ -12,10 +47,11 @@ %% How to input data? An example. (cf. test13.sm1) %% Modify the following or copy the /cone.sample { ... } def %% to your own file, -%% edit it, and execute if by " cone.sample ; " +%% edit it, and execute it by " cone.sample ; " %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% /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 @@ -75,6 +111,7 @@ def ] def +/cone.DhH 0 def % Set a function to compute Grobner basis. % cone.gb_Dh : For computing in Homogenized Weyl algebra h[1,1](D). % cone.gb_DhH : For computing in doubly homogenized Weyl algebra. @@ -88,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 @@ -96,7 +136,7 @@ getGrobnerFan printGrobnerFan %%%% If you want to save the data to the file sm1out.txt, then uncomment. -% /cone.wightGblist 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). @@ -127,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 @@ -193,6 +234,7 @@ def cone.input { . homogenize toString } map /cone.input set dh.end +/cone.DhH 1 def % Set a function to compute Grobner basis. % cone.gb_Dh : For computing in Homogenized Weyl algebra h[1,1](D). % cone.gb_DhH : For computing in doubly homogenized Weyl algebra. @@ -213,7 +255,7 @@ getGrobnerFan printGrobnerFan %%%% If you want to save the data to the file sm1out.txt, then uncomment. -% /cone.wightGblist 1 def saveGrobnerFan /ff set ff output +% /cone.withGblist 1 def saveGrobnerFan /ff set ff output %%%% Step 2. Dehomogenize the Grobner Cones %%%% by the equivalence relation in a local ring (uncomment). @@ -243,18 +285,18 @@ dhcone.printGrobnerFan % If you use local polymake, then comment out. % If you use the cgi/polymake on the net, then uncomment out. -%/doPolymake {doPolymake.OoHG} def - -% My setting. -[(getenv) (HOST)] extension /cone.hostname set -cone.hostname tag 0 eq { /cone.hostname (?) def } { } ifelse -cone.hostname (orange2.math.sci.kobe-u.ac.jp) eq -cone.hostname (orange2-clone.math.sci.kobe-u.ac.jp) eq +%/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.web 1 eq or -{ - (Using doPolymake.OoHG ) message - /doPolymake {doPolymake.OoHG} def -} { } ifelse +{ + (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 @@ -267,10 +309,12 @@ or /cone.loaded boundp { } { [(parse) (cohom.sm1) pushfile] extension -% [(parse) (cone.sm1) pushfile] extension +% [(parse) (cone.sm1) pushfile] extension % BUG? cone.sm1 overrides a global + % in cohom.sm1? [(parse) (dhecart.sm1) pushfile] extension /cone.loaded 1 def - oxNoX polymake.start ( ) message + oxNoX + polymake.start ( ) message } ifelse } def @@ -290,8 +334,10 @@ or /dh.autoHomogenize 0 def [(AutoReduce) 1] system_variable [ff { toString } map cone.vv - [ww cone.vv generateD1_1]] dh.gb 0 get /arg1 set + [ww cone.vv generateD1_1]] ff getAttributeList setAttributeList + dh.gb 0 get /arg1 set ] pop + popVariables arg1 } def @@ -395,8 +441,36 @@ cone.comment message % cone.d (pointed cones lies in this space. cf. cone.Lp) % These are set during getting the cone.startingCone +%< +% global +%cone.ckmFlip. Collar-Kalkbrener-Mall の flip アルゴリズムを使わない 0. 使う 1. +% Default は 0. +%> +/cone.ckmFlip 0 def %< +% global +% cone.DhH dx x = x dx + h H なら 1. dx x = x dx + h^2 なら 0. Default 0. +%> +/cone.DhH 0 def + +%< +% Global +% gbCheck をするか? しないと結果はあやふや. しかしメモリ exhaust は防げる. +% 使うときは /cone.epsilon, /cone.epsilon.limit を十分小さくしておく. +%> +/cone.do_gbCheck 1 def + +% Default の cone.gb の定義. 各プログラムで再度定義してもよい. +/cone.gb { + cone.DhH { + cone.gb_DhH + } { + cone.gb_Dh + } ifelse +} def + +%< % Usage: wv g coneEq1 % in(f) が monomial 専用. in_w(f) = LT(f) となる weight w の満たす % 不等式制約を求める. @@ -985,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 @@ -1010,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 @@ -1586,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 [ @@ -1843,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 @@ -2020,14 +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] 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 @@ -2323,22 +2433,25 @@ def } def %< -% usages: getNextFlip [cone, k] +% usages: getNextFlip [cone, k, cid] % cone.fan を検索して まだ flip してない cone と facet の組を戻す. % もうないときには null を戻す. +% cid は cone が cone.fan の 何番目であるかの index. cone.gblist の検索等に +% 用いる. %> /getNextFlip { - [/tcone /ans /ii ] pushVariables + [/tcone /ans /ii /cid] pushVariables [ - /ans null def + /ans null def /cid -1 def 0 1 cone.fan length 1 sub { /ii set cone.fan ii get /tcone set + /cid ii def tcone getNextFacet /ans set ans tag 0 eq { } { exit } ifelse } for ans tag 0 eq { /arg1 null def } - { /arg1 [tcone ans] def } ifelse + { /arg1 [tcone ans cid] def } ifelse ] pop popVariables arg1 @@ -2348,12 +2461,13 @@ def % flip の時の epsilon /cone.epsilon (1).. (10).. div def /cone.epsilon.limit (1).. (100).. div def +% cone.epsilon.limit を負にすれば停止しない. %< % Usages: result_getNextFlip getNextCone ncone % flip して新しい ncone を得る. %> -/getNextCone { +/getNextCone.orig { /arg1 set [/ncone /ccone /kk /w /next_weight_w_wv] pushVariables [ @@ -2687,6 +2801,7 @@ def /printGrobnerFan { [/i] pushVariables [ + $(gfan) usage to find explanations on variables.$ message (========== Grobner Fan ====================) message [ (cone.comment) @@ -3412,3 +3527,407 @@ def % % Todo: save functions. + +%< +% Collart, Kalkbrener, Mall のアルゴリズムによる gb の flip. +% See also Sturmfels' book, p.22, 23. +% Usages: [reducedGb, vlist, oldWeight, facetWeight, newWeight] ckmFlip rGb +% If it fails, then it returns null, else it returns the reducedGb for the +% newWeight. +% gb の check をやるので, それに失敗したら null を戻す. +% weight はすべて vw 形式で. vw 形式 = variable weight の繰り返しの形式 +% reducedGb は文字列のリストではなく多項式の形式のこと. +% 理由は reducedGb より ring の構造を読むため. +%> +/ckmFlip { + /arg1 set + [/arg_ckmFlip /gOld /vlist /oldWeight /facetWeight /newWeight + /gNew + /ww /ww1 /ww2 % 本の中の w1, w, w2 (古い, facet, 新しい) + /ch1 /ch2 % 本の中の {\cal H}_1, {\cal H}_2 + /grData /rTable + /rTable2 % rTable の反対の変換. + /facetWeight_gr /vlist_gr % graded ring 用. + /oldWeight_gr + /ccf % reduction した係数. + /rwork /ccf2 /gNew + ] pushVariables + [ + arg1 /arg_ckmFlip set + arg_ckmFlip 0 get /gOld set + arg_ckmFlip 1 get /vlist set + arg_ckmFlip 2 get /oldWeight set + arg_ckmFlip 3 get /facetWeight set + arg_ckmFlip 4 get /newWeight set + +% facet weight vector ww についての initial を取り出す. ch1 へいれる. + gOld getRing ring_def + facetWeight weightv /ww set + gOld { ww init } map /ch1 set % facetWeight による initial の取り出し. + + +% 例: [(x,y) [(x) -1 (Dx) 1 (y) -1 (Dy) 2]] getGrRing +% [$x,y,y',$ , [ $x$ , $y$ ] , [ [ $Dy$ , $y'$ ] ] ] +% 変数リスト 置換表 +% ch1 を gr_ww の元に変換. + [vlist facetWeight] getGrRing /grData set + [grData 0 get ring_of_differential_operators 0] define_ring /rwork set + grData 2 get { { . } map } map /rTable set + rTable { reverse } map /rTable2 set + grData 0 get /vlist_gr set + ch1 { toString . rTable replace toString } map /ch1 set + + oldWeight { dup isString { . rTable replace toString } + { } ifelse } map /oldWeight_gr set + +% facetWeight も 新しい環 gr_ww の weight に変換. +% 例. [(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 + %% ttt cone.gb_gr /ch1 set %再度の計算は不要. + [[(1)] vlist_gr oldWeight_gr] cone.gb_gr getRing ring_def % Set Ring. + ch1 {toString .} map /ch1 set +%% ここまででとりあえずテストをしよう. +%% ch1 /arg1 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 で計算. +% どちらをとるかは cone.reduction_gr で区別するしかなし + ch1 getRing ring_def ; + ch2 {toString .} map {ch1 cone.reduction} map /ccf set + %ccf pmat + % とりあえずテスト. + % [ch1 ch2] /arg1 set + %% ccf[i][0] は 0 でないと矛盾. check まだしてない. + + %% ccf[i][2] (syzygy) を gr から もとの ring へ戻し, + %% 新しい reduced gbasis を ccf[i][2] * gOld で作る. + rwork ring_def + ccf { 2 get {toString . rTable2 replace toString} map } map /ccf2 set + %% ccf2 は gr でない ring の元. + gOld getRing ring_def + cone.DhH { cone.begin_DhH } { } ifelse % Hh か h^2 か. + ccf2 { {.} map gOld mul } map /gNew set + gNew { toString } map /gNew set + cone.DhH { cone.end_DhH } { } ifelse % Hh か h^2 か. + % gNew /arg1 set + %gNew が newWeight での GB か check. Yes なら reduced basis へ. + %No なら null を戻す. +%%Ref: note @s/2005/06/30-note-gfan.pdf + 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 + ] pop + popVariables + arg1 +} def + +%< +% Usages: f gbasis cone.reduction_DhH +% dx x = x dx + h H での reduction. +%> +/cone.reduction_DhH { + /arg2 set /arg1 set + [/ff /ggbasis /eenv /ans] pushVariables + [ + /ff arg1 def /ggbasis arg2 def + cone.begin_DhH + ff ggbasis reduction /ans set + cone.end_DhH + /arg1 ans def + ] pop + popVariables + arg1 +} def + +%< +% Usages: f gbasis cone.reduction_Dh +% dx x = x dx + h^2 での reduction. +%> +/cone.reduction_Dh { + /arg2 set /arg1 set + [/ff /ggbasis /eenv /ans] pushVariables + [ + /ff arg1 def /ggbasis arg2 def + ff ggbasis reduction /ans set + /arg1 ans def + ] pop + popVariables + arg1 +} def + +%< +% Usages: cone.begin_DhH dx x = x dx + h H を開始. +%> +/cone.begin_DhH { + [(Homogenize) (AutoReduce) (KanGBmessage)] pushEnv /cone.eenv set + [(Homogenize) 3] system_variable +} def + +%< +% Usages: cone.begin_DhH dx x = x dx + h H を終了. +%> +/cone.end_DhH { + cone.eenv popEnv +} def + +%< +% Usages: ff vv ww cone.gb_gr_DhH dx x = x dx + h H で計算. +% dh.gb は dhecart.sm1 で定義されており, dx x = x dx + h H での計算. +% gr をとっても, -w,w の場合は 微分作用素環のままであり, これが必要. +% bug? cone.gb で十分? +%> +/cone.gb_gr_DhH { + /arg1 set + [/ff /ww /vv] pushVariables + [ + /ff arg1 0 get def + /vv arg1 1 get def + /ww arg1 2 get def + /dh.gb.verbose 1 def + /dh.autoHomogenize 0 def + [(AutoReduce) 1] system_variable + [ff { toString } map vv + [ww vv generateD1_1]] dh.gb 0 get /arg1 set + ] pop + popVariables + arg1 +} def +%< +% Usages: ff vv ww cone.gb_gr_Dh dx x = x dx + h^2 で計算. +% gb は dhecart.sm1 で定義されており, dx x = x dx + h^2 での計算. +% gr をとっても, -w,w の場合は 微分作用素環のままであり, これが必要. +% bug? cone.gb で十分? +%> +/cone.gb_gr_Dh { + /arg1 set + [/ff /ww /vv /gg /envtmp] pushVariables + [ + /ff arg1 0 get def + /vv arg1 1 get def + /ww arg1 2 get def + + [(AutoReduce) (KanGBmessage)] pushEnv /envtmp set + [(AutoReduce) 1] system_variable + [(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 +} def + + +% これらは cone.ckmFlip 1 の時しか使わず. +/cone.reduction { + cone.DhH { + cone.reduction_DhH + }{ + cone.reduction_Dh + } ifelse +} def +/cone.gb_gr { + cone.DhH { + cone.gb_gr_DhH + }{ + cone.gb_gr_Dh + } ifelse +} def + + +/test1.ckmFlip { + % cf. cone.sample2 + cone.load.cohom + /cone.comment [ + (BS for y and y-(x-1)^2, t1, t2 space, in doubly homogenized Weyl algebra.) nl + (The Grobner cones are dehomogenized to get local Grobner fan.) nl + ] cat def + /cone.vlist [(t1) (t2) (x) (y) (Dt1) (Dt2) (Dx) (Dy) (h) (H)] def + /cone.vv (t1,t2,x,y) def + /cone.type 1 def + /cone.parametrizeWeightSpace { + 4 2 parametrizeSmallFan + } def + + /cone.DhH 1 def + /cone.ckmFlip 1 def + + /cone.local 1 def + /cone.w_start null def + /cone.h0 1 def + /cone.input + [ + (t1-y) (t2 - (y-(x-1)^2)) + ((-2 x + 2)*Dt2+Dx) + (Dt1+Dt2+Dy) + ] + def + % homogenize + [cone.vv ring_of_differential_operators + [[(t1) -1 (t2) -1 (Dt1) 1 (Dt2) 1]] ecart.weight_vector + 0] define_ring + dh.begin + cone.input { . homogenize toString } map /cone.input set + dh.end + + +% テストを開始する. +% getStartingCone /cone.ncone set +% cone.ncone updateFan +% cone.gblist 0 get message +% cone.ncone /cone.ccone set +% getNextFlip /cone.nextflip set +% cone.nextflip message + + /wOld [(t1) , -29 , (t2) , -38 , (Dt1) , 29 , (Dt2) , 38 ] def + /wFacet [(t1) , -1 , (t2) , -1 , (Dt1) , 1 , (Dt2) , 1 ] def + /wNew [(t1) , -39 , (t2) , -38 , (Dt1) , 39 , (Dt2) , 38 ] def + cone.input wOld cone.gb /ff set + [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