[BACK]Return to gfan.sm1 CVS log [TXT][DIR] Up to [local] / OpenXM / src / kan96xx / Doc

File: [local] / OpenXM / src / kan96xx / Doc / gfan.sm1 (download)

Revision 1.4, Wed Sep 15 07:41:59 2004 UTC (19 years, 9 months ago) by takayama
Branch: MAIN
Changes since 1.3: +366 -6 lines

New function  dhCones_h  dehomogenizes cones (cone.fan) obtained
by getGrobnerFan in the doubly homogenized ring (h, H).
Cones which have same initials in the (local) graded ring are merged into
one ring.
The result is stored in the variable dhcone.fan.

%  $OpenXM: OpenXM/src/kan96xx/Doc/gfan.sm1,v 1.4 2004/09/15 07:41:59 takayama Exp $
% cp cone.sm1 $OpenXM_HOME/src/kan96xx/Doc/gfan.sm1
% $Id: cone.sm1,v 1.48 2004/09/15 07:38:42 taka Exp $
% iso-2022-jp

/cone.debug 1 def

/ox.k0.loaded boundp {
} { 
 [(parse) (ox.sm1) pushfile] extension
} ifelse

%
% cone.fan, cone.gblist $B$K(B fan $B$N%G!<%?$,$O$$$k(B.
%

%%%%<<<<  $B=i4|%G!<%?$N@_DjNc(B  data/test13 $B$h$j(B.  <<<<<<<<<<<<<<
/cone.sample.test13 {
 /cone.loaded boundp { }
 {
  [(parse) (cohom.sm1) pushfile] extension
  [(parse) (cone.sm1) pushfile] extension
  /cone.loaded 1 def
 } ifelse
/cone.comment [
  (Toric ideal for 1-simplex x 2-simplex, in k[x]) nl
] cat def
%------------------Globals----------------------------------------
% Global: cone.type
% $B$I$N(B exponents $B$r<h$j=P$9$N$+;XDj$9$k(B.
% cf. exponents, gbext  h $B$d(B H $B$b8+$k$+(B?
% 0 : x,y,Dx,Dy
% 1 : x,y,Dx,Dy,h,H
% 2 : x,y,Dx,Dy,h 
/cone.type 2 def

% Global: cone.local
% cone.local: Local $B$+(B?  1 $B$J$i(B local
/cone.local 0 def	


% Global: cone.h0
% cone.h0:  1 $B$J$i(B h $B$N(B weight 0 $B$G$N(B Grobner fan $B$r7W;;$9$k(B.
/cone.h0 1 def

% ---------------  $BF~NO%G!<%?MQBg0hJQ?t$N@_Dj(B --------------------------
%
% cone.input : $BF~NOB?9`<07O(B
/cone.input 
  [
    (x11 x22 - x12 x21) (x12 x23 - x13 x22)
    (x11 x23 - x13 x21)
  ]
def

% cone.vlist : $BA4JQ?t$N%j%9%H(B
/cone.vlist [(x11) (x12) (x13) (x21) (x22) (x23)
             (Dx11) (Dx12) (Dx13) (Dx21) (Dx22) (Dx23) (h)] def

% cone.vv : define_ring $B7A<0$NJQ?t%j%9%H(B.
/cone.vv (x11,x12,x13,x21,x22,x23) def

% cone.parametrizeWeightSpace : weight $B6u4V$r(B parametrize $B$9$k4X?t(B. 
%   $BBg0hJQ?t(B cone.W , cone.Wpos $B$b$-$^$k(B.
/cone.parametrizeWeightSpace {
  6 6 parametrizeSmallFan
} def

% cone.w_start : weight$B6u4V$K$*$1$k(B weight $B$N=i4|CM(B. 
% $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. 
% random $B$K$d$k$H$-$O(B null $B$K$7$F$*$/(B.
/cone.w_start
  [9 8 5 4 5 6]
def   

% cone.gb : gb $B$r7W;;$9$k4X?t(B.
/cone.gb {
  cone.gb_Dh
} def



( ) message
cone.comment message
(cone.input = ) messagen cone.input message
(Type in getGrobnerFan) message
(Do clearGlobals if necessary) message
(printGrobnerFan ; saveGrobnerFan /ff set ff output ) message

} def
%%%%%%>>>>>  $B=i4|%G!<%?$N@_DjNc$*$o$j(B >>>>>>>>>>>>>>>>>>>>>>

% Global: cone.type
% $B$I$N(B exponents $B$r<h$j=P$9$N$+;XDj$9$k(B.
% cf. exponents, gbext  h $B$d(B H $B$b8+$k$+(B?
% 0 : x,y,Dx,Dy
% 1 : x,y,Dx,Dy,h,H
% 2 : x,y,Dx,Dy,h 
/cone.type 2 def

% Global: cone.local
% cone.local: Local $B$+(B?  1 $B$J$i(B local
/cone.local 1 def

% Global: cone.h0
% cone.h0:  1 $B$J$i(B h $B$N(B weight 0 $B$G$N(B Grobner fan $B$r7W;;$9$k(B.
/cone.h0 1 def

% Global: cone.n (number of variables in GB)
%         cone.m (freedom of the weight space. cf. cone.W)
%         cone.d (pointed cones lies in this space. cf. cone.Lp)
% These are set during getting the cone.startingCone


%<
% Usage:  wv g coneEq1
% in(f) $B$,(B monomial $B@lMQ(B.  in_w(f) = LT(f) $B$H$J$k(B weight w $B$NK~$?$9(B
% $BITEy<0@)Ls$r5a$a$k(B.
%>
/coneEq1 {
  /arg1 set
  [/g /eqs /gsize /i /j /n /f /exps /m  % Do not use "eq" as a variable
   /expsTop
  ] pushVariables
  [
    /g arg1 def  % Reduced Grobner basis
    /eqs [ ] def % $BITEy<07O$N78?t(B
    /gsize g length def
    0 1 gsize 1 sub {
      /i set
      g i get /f set  % f $B$O(B i $BHVL\$N(B reduced Grobner basis $B$N85(B
      [(exponents) f cone.type] gbext /exps set % exps $B$O(B f $B$N(B exponent vector
      exps length /m set
      m 1 eq not {
        /expsTop exps 0 get def % expsTop $B$O(B f $B$N@hF,$N(B exponent vector.
        1 1 exps length 1 sub {
           /j set
           eqs [expsTop exps j get  sub] join /eqs set 
           % exps[0]-exps[j] $B$r(B eqs $B$X3JG<$7$F$$$/$@$1(B.
           % Cone $B$N(B closure $B$r$@$9$N$G(B  >= $B$G(B OK.
        } for
      } { } ifelse
    } for
    /arg1 eqs def   
  ] pop
  popVariables 
  arg1
} def

%<
% Usage: ww g coneEq
% 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. 
% g $B$O(B reduced Grobner basis
% in(f) $B$,(B monomial $B$G$J$$>l9g$b07$&(B.
% in_w(f) = in_ww(f) $B$H$J$k(B weight w $B$NK~$?$9(B
% $BITEy<0@)Ls$r5a$a$k(B.
% ord_w, init (weightv) $B$rMQ$$$k(B.
%>
/coneEq {
  /arg2 set
  /arg1 set
  [/g /eqs /gsize /i /j /n /f /exps /m 
   /expsTop /ww /ww2 /iterms
  ] pushVariables
  [
    /g arg2 def  % Reduced Grobner basis
    /ww arg1 def % weight vector. v-w $B7A<0(B
    ww to_int32 /ww set % univNum $B$,$"$l$P(B int32 $B$KD>$7$F$*$/(B.
    /ww2 ww weightv def  % v-w $B7A<0$r(B $B?t;z$N%Y%/%H%k$K(B. (init $BMQ(B)

    /eqs null def % $BITEy<07O$N78?t(B
    /gsize g length def
    0 1 gsize 1 sub {
      /i set
      g i get /f set  % f $B$O(B i $BHVL\$N(B reduced Grobner basis $B$N85(B
      [(exponents) f cone.type] gbext /exps set % exps $B$O(B f $B$N(B exponent vector
      exps length /m set
      m 1 eq not {
        /expsTop exps 0 get def % expsTop $B$O(B f $B$N@hF,$N(B exponent vector.
        /iterms f ww2 init length def % f $B$N(B initial term $B$N9`$N?t(B.
        % in_ww(f) > f_j $B$H$J$k9`$N=hM}(B.
        iterms 1 exps length 1 sub {
           /j set
           expsTop exps j get sub    eqs cons /eqs set 
           % exps[0]-exps[j] $B$r(B eqs $B$X3JG<$7$F$$$/(B.
        } for
        % in_ww(f) = f_j $B$H$J$k9`$N=hM}(B.
        [(exponents) f ww2 init cone.type] gbext /exps set % exps $B$O(B in(f)
        1 1 iterms 1 sub {
          /j set
          exps j get expsTop sub   eqs cons /eqs set 
          expsTop exps j get sub   eqs cons /eqs set 
          % exps[j]-exps[0], exps[0]-exps[j] $B$r3JG<(B. 
          % $B7k2LE*$K(B (exps[j]-exps[0]).w = 0 $B$H$J$k(B.  
        }  for
      } { } ifelse
    } for
    eqs listToArray reverse /eqs set
    /arg1 eqs def   
  ] pop
  popVariables 
  arg1
} def

%<
% Usage: wv g coneEq genPo
% polymake $B7A<0$N(B INEQUALITIES $B$r@8@.$9$k(B.  coneEq -> genPo $B$HMxMQ(B
%>
/genPo {
  /arg1 set
  [/outConeEq /rr /nn /ii /mm /jj /ee] pushVariables
  [
    /outConeEq arg1 def
    /rr [(INEQUALITIES) nl] cat def % $BJ8;zNs(B rr $B$KB-$7$F$$$/(B.
    outConeEq length /nn set
    0 1 nn 1 sub {
      /ii set
      outConeEq ii get /ee set
      [ rr 
        (0 )    % $BHs$;$$$8MQ$N(B 0 $B$r2C$($k(B.
        0 1 ee length 1 sub {
          /jj set
          ee jj get toString ( )
        } for
        nl 
      ] cat /rr set       
    } for
    /arg1 rr def
  ] pop
  popVariables
  arg1
} def

%<
% Usage: wv g coneEq genPo2
% doPolyamke $B7A<0$N(B INEQUALITIES $B$r@8@.$9$k(B.  coneEq -> genPo2 $B$HMxMQ(B
% tfb $B7A<0J8;zNs(B.
%>
/genPo2 {
  /arg1 set
  [/outConeEq /rr /nn /ii /mm /jj /ee] pushVariables
  [
    /outConeEq arg1 def
    /rr $polymake.data(polymake.INEQUALITIES([$ def 
    % $BJ8;zNs(B rr $B$KB-$7$F$$$/(B.
    outConeEq length /nn set
    0 1 nn 1 sub {
      /ii set
      outConeEq ii get /ee set
      [ rr 
        ([0,)   % $BHs$;$$$8MQ$N(B 0 $B$r2C$($k(B.
        0 1 ee length 1 sub {
          /jj set
          ee jj get toString 
          jj ee length 1 sub eq { } { (,) } ifelse
        } for
        (]) 
        ii nn 1 sub eq { } { (,) } ifelse
      ] cat /rr set       
    } for
    [rr $]))$ ] cat /rr set
    /arg1 rr def
  ] pop
  popVariables
  arg1
} def

/test1 {
  [(x,y) ring_of_differential_operators 0] define_ring
  [ (x + y + Dx + Dy).
    (x ^2 Dx^2 + y^2 Dy^2).
    (x). 
  ] /gg set
  gg coneEq1 /ggc set
  gg message
  ggc pmat 

  ggc genPo message
} def

/test2 {
  [(parse) (dhecart.sm1) pushfile] extension
  dh.test.p1 /ff set 
  ff 0 get coneEq1 /ggc set
  ggc message
  ggc genPo /ss set
  ss message
  (Data is in ss) message
} def


/test3 {
%  [(parse) (cohom.sm1) pushfile] extension
  /ww [(Dx) 1 (Dy) 1] def
  [(x,y) ring_of_differential_operators 
   [ww] weight_vector
   0] define_ring
  [ (x Dx + y Dy -1).
    (y^2 Dy^2 + 2 + y Dy ).
  ] /gg set
  gg {homogenize} map /gg set
  [gg] groebner 0 get /gg set
  ww message 
  ww gg coneEq /ggc set
  gg message
  ggc pmat 

  ggc genPo message
} def

%<
% Usage: test3b
% Grobner cone $B$r7hDj$7$F(B, polymake $BMQ$N%G!<%?$r@8@.$9$k%F%9%H(B.
% weight (0,0,1,1) $B$@$H(B max dim cone $B$G$J$$(B.
%>
/test3b {
%  [(parse) (cohom.sm1) pushfile] extension
  /ww [(Dx) 1 (Dy) 2] def
  [(x,y) ring_of_differential_operators 
   [ww] weight_vector
   0] define_ring
  [ (x Dx + y Dy -1).
    (y^2 Dy^2 + 2 + y Dy ).
  ] /gg set
  gg {homogenize} map /gg set
  [gg] groebner 0 get /gg set
  ww message 
  ww gg coneEq /ggc set
  gg message
  ggc pmat 

%  ggc genPo /ggs set % INEQ $B$rJ8;zNs7A<0$G(B
%  ggs message
%  ggs output
%  (mv sm1out.txt test3b.poly) system
%  (Type in polymake-pear.sh test3b.poly FACETS) message

   ggc genPo2 /ggs set % INEQ $B$rJ8;zNs7A<0(B for doPolymake
   ggs message

} def

% commit (dr.sm1):  lcm, denominator, ngcd, to_univNum,  numerator, reduce 
%  8/22,  changelog-ja $B$^$@(B. 
% to do : nnormalize_vec,  sort_vec --> shell $B$G(B OK.
% 8/27, getNode

/test4 {
 $polymake.data(polymake.INEQUALITIES([[0,1,0,0],[0,0,1,0]]))$ /ff set 
 [(FACETS) ff] doPolymake /rr set

 rr 1 get /rr1 set
 rr1 getLinearitySubspace pmat

} def

%<
% Usage: vv ineq isInLinearSpace
%        vv $B$,(B ineq[i] > 0 $B$GDj5A$5$l$kH>6u4V$N$I$l$+$K$O$$$C$F$$$k$J$i(B 0
%        vv $B$,(B $BA4$F$N(B i $B$K$D$$$F(B ineq[i] = 0 $B$K$O$$$C$F$$$?$i(B 1.
%>
/isInLinearSpace {
  /arg2 set
  /arg1 set
  [/vv /ineq /ii /rr] pushVariables
  [
    /vv arg1 def
    /ineq arg2 def
    /rr 1 def
    {
       0 1 ineq length 1 sub {
         /ii set
         % vv . ineq[ii] != 0 $B$J$i(B vv $B$O(B linearity space $B$N85$G$J$$(B.
         vv ineq ii get mul to_univNum isZero {
         } { /rr 0 def exit} ifelse
       } for
       exit
    } loop
    /arg1 rr def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: doPolymakeObj getLinearitySubspace
% INEQUALITIES $B$H(B VERTICES $B$+$i(B maximal linearity subspace
% $B$N@8@.%Y%/%H%k$r5a$a$k(B.
% $BNc(B: VERTICES [[0,1,0,0],[0,0,1,0],[0,0,0,-1],[0,0,0,1]]] 
% $BNc(B: INEQUALITIES [[0,1,0,0],[0,0,1,0]]
% $BF~NO$O(B polymake $B$N(B tree (doPolymake $B$N(B 1 get)
%>
/getLinearitySubspace {
  /arg1 set
  [/pdata /vv /ineq /rr /ii] pushVariables
  [
     /pdata arg1 def
     {
       /rr [ ] def
       % POINTED $B$J$i(B max lin subspace $B$O(B 0.
       pdata (POINTED) getNode tag 0 eq { } { exit} ifelse 

       pdata (INEQUALITIES) getNode 2 get 0 get /ineq set
       pdata (VERTICES) getNode 2 get 0 get /vv set
       0 1 vv length 1 sub {
         /ii set
         % -vv[ii] $B$,(B ineq $B$rK~$?$9$+D4$Y$k(B. 
         vv ii get ineq  isInLinearSpace {
            rr  [vv ii get] join /rr set 
         } {  } ifelse
       } for
       exit
     } loop
     /arg1 rr def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: mm asir_matrix_image
% $B@8@.85$h$j@~7A6u4V$N4pDl$rF@$k(B.
%>
/asir_matrix_image {
  /arg1 set
  [/mm /rr] pushVariables
  [(CurrentRingp)] pushEnv
  [
    /mm arg1 def
    mm to_univNum /mm set
    oxasir.ccc [ ] eq {
       (Starting ox_asir server.) message
        ox_asirConnectMethod
    } {  } ifelse
    {
     oxasir.ccc [(matrix_image) mm] asir
     /rr set
     rr null_to_zero /rr set
     exit

     (asir_matrix_image: not implemented) error exit
    } loop

    rr numerator /rr set
    /arg1 rr def
  ] pop
  popEnv
  popVariables
  arg1
} def
[(asir_matrix_image)
 [(Calling the function matrix_image of asir. It gets a reduced basis of a given matrix.)
  (Example:  [[1 2 3] [2 4 6]] asir_matrix_image)
]] putUsages

%<
% Usages: mm asir_matrix_kernel
% $BD>8r$9$k6u4V$N4pDl(B.
%>
/asir_matrix_kernel {
  /arg1 set
  [/mm /rr] pushVariables
  [(CurrentRingp)] pushEnv
  [
    /mm arg1 def
    mm to_univNum /mm set
    oxasir.ccc [ ] eq {
       (Starting ox_asir server.) message
        ox_asirConnectMethod
    } {  } ifelse
    {
     oxasir.ccc [(matrix_kernel) mm] asir
     /rr set
     rr null_to_zero /rr set
     exit

     (asir_matrix_image: not implemented) error exit
    } loop
    rr 1 get numerator /rr set  
    /arg1 rr def
  ] pop
  popEnv
  popVariables
  arg1
} def
[(asir_matrix_kernel)
 [(Calling the function matrix_kernel of asir.)
  (It gets a reduced basis of the kernel of a given matrix.)
  (Example:  [[1 2 3] [2 4 6]] asir_matrix_kernel)
]] putUsages

%<
% Usages: v null_to_zero
%>
/null_to_zero {
  /arg1 set
  [/pp /rr] pushVariables
  [
    /pp arg1 def
    {
      /rr pp def
      pp isArray {
       pp {null_to_zero} map /rr set
       exit
      }{ } ifelse
    
      pp tag 0 eq {
        /rr (0).. def
        exit
      }{  } ifelse
      exit
    } loop
    /arg1 rr def
  ] pop
  popVariables
  arg1 
} def
[(null_to_zero)
[(obj null_to_zero rob)
 $It translates null to (0)..$
]] putUsages

%<
% Usages: newVector.with-1  
% (-1).. $B$GKd$a$?%Y%/%H%k$r:n$k(B.
%>
/newVector.with-1 {
  newVector { pop (-1).. } map
} def


% [2 0] lcm $B$O(B 0 $B$r$b$I$9$,$$$$$+(B? --> OK.

%<
% Usages: mm addZeroForPolymake
% $B0J2<$NFs$D$N4X?t$O(B,  toQuotientSpace $B$K$bMxMQ(B.
% Polymake INEQUALITIES $BMQ$K(B 0 $B$r;O$a$KB-$9(B.
% $BF~NO$O(B $B%j%9%H$N%j%9%H(B
% [[1,2], [3,4],[5,6]] --> [[0,1,2],[0,3,4],[0,5,6]]
%>
/addZeroForPolymake {
  /arg1 set
  [/mm /rr] pushVariables
  [
    /mm arg1 def
    mm to_univNum /mm set
    mm { [(0)..] 2 1 roll join } map /mm set
    /arg1 mm def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: mm cone.appendZero
%>
/cone.appendZero {
  /arg1 set
  [/mm /rr] pushVariables
  [
    /mm arg1 def
    mm to_univNum /mm set
    mm { [(0)..] join } map /mm set
    /arg1 mm def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: mm removeFirstFromPolymake
% $B;O$a$N(B 0 $B$r<h$j=|$/(B.
% $BF~NO$O(B $B%j%9%H$N%j%9%H(B
% [[0,1,2],[0,3,4],[0,5,6]] ---> [[1,2], [3,4],[5,6]]
%>
/removeFirstFromPolymake {
  /arg1 set
  [/mm /rr] pushVariables
  [
    /mm arg1 def
    mm to_univNum /mm set
    mm {rest} map /mm set
    /arg1 mm def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: mm genUnit
% [1,0,0,...] $B$r2C$($k$?$a$K@8@.(B.
% [[0,1,2], [0,3,4],[0,5,6]]--> [1,0,0]
%>
/genUnit {
  /arg1 set
  [/mm /rr /i] pushVariables
  [
    /mm arg1 def
    mm 0 get length newVector /rr set
    rr null_to_zero /rr set
    rr 0 (1).. put
    /arg1 rr def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: mm genUnitMatrix
% [[0,1,2], [0,3,4],[0,5,6]]--> [[1,0,0],[0,1,0],[0,0,1]]
%>
/genUnitMatrix {
  /arg1 set
  [/mm /rr /nn /i] pushVariables
  [
    /mm arg1 def
    mm 0 get length /nn set
    [
      0 1 nn 1 sub {
        /i set
        nn newVector null_to_zero /mm set 
        mm i (1).. put
        mm
      } for
    ]
    /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
%%note:  2004, 8/29 (sun)
% toQuotientSpace : Linearity space $B$G3d$k(B.
% Usages: ineq mm toQuotientSpace
% $BF~NO$O(B coneEq $B$N=PNO(B ineq  
% $B$*$h$S(B doPolymake --> getLinearitySubspace ==> L 
%  [L,[1,0,0,...]] asir_matrix_kernel removeFirstFromPolymake $B$GF@$i$l$?(B mm
% $B=PNO$+$i(B 0 $B%Y%/%H%k$O:o=|(B.
% $B=PNO$b(B coneEq $B7A<0(B.  $BFC$K(B polymake $BMQ$K(B 0 $B$r2C$($k$N$,I,MW(B.
% ref: getUnit, removeFirstFromPolymake, addZeroForPolymake,
%      asir_matrix_kernel, getLinearitySubspace
%>
/toQuotientSpace {
  /arg2 set
  /arg1 set
  [/ineq /mm /rr] pushVariables
  [
    /ineq arg1 def
    /mm arg2 def

    ineq mm transpose mul /rr set

    /arg1 rr def
  ] pop
  popVariables
  arg1
} def

/test5.data
 $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]]))$
def
%<
% Usages: test5
%% getConeInfo $B$rJQ99$9$l$P(B polymake $B$r8F$P$:$K%F%9%H$G$-$k(B.
%>
/test5 {
  % test3b $B$h$j(B
  /ww [(Dx) 1 (Dy) 2] def
%  /ww [(x) 1 (y) -2 (Dx) 3 (Dy) 6] def
  [(x,y) ring_of_differential_operators 
   [ww] weight_vector
   0] define_ring
  [ (x Dx + y Dy -1).
    (y^2 Dy^2 + 2 + y Dy ).
  ] /gg set
  gg {homogenize} map /gg set
  [(AutoReduce) 1] system_variable
  [gg] groebner 0 get /gg set
  ww message 

  ww gg coneEq getConeInfo /rr set
  (Type in rr 0 get :: ) message
} def
%[5, [[1,0,1,0,-2],[0,1,0,1,-2]], $NOT__POINTED$ ]
% $B$3$N>l9g$O(B 2 $B<!85$^$GMn$9$H(B pointed cone $B$K$J$k(B.
%  coneEq mmc transpose $B$r$b$H$K(B FACETS $B$r7W;;$9$l$P$h$$(B.

%<
% Usage: ceq getConeInfo
% 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. 
% 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.
% Grobner cone $B$N(B $B<!85(B cdim (DIM), $BJd6u4V(B (linearity space ) $B$X$N9TNs(B mmc
% linearity space $B<+BN(B, pointed or not__pointed
% $B$D$^$j(B [cdim, L', L, PointedQ]
% $B$r7W;;$7$FLa$9(B.  (polymake $B7A<0$NM>J,$JItJ,$J$7(B)
% polymake $BI,MW(B.
% ref: coneEq
% Global:
% cone.getConeInfo.rr0, cone.getConeInfo.rr1 $B$K(B polymake $B$h$j$NLa$jCM$,$O$$$k(B.
%>
/getConeInfo {
  /arg1 set
  [/ww /g /ceq /ceq2 /cdim /mmc /mmL /rr /ineq /ppt] pushVariables
  [
     /ceq arg1 def 
     ceq pruneZeroVector /ceq set
     ceq genPo2 /ceq2 set
     % ceq2 $B$O(B polymake.data(polymake.INEQUALITIES(...)) $B7A<0(B
     % polymake $B$G(B ceq2 $B$N<!85$N7W;;(B.
     /getConeInfo.ceq  ceq def /getConeInfo.ceq2 ceq2 def

     cone.debug { (Calling polymake DIM.) message } { } ifelse
     [(DIM) ceq2] doPolymake 1 get /rr set
     cone.debug {(Done.) message } {  } ifelse
% test5 $B$K$O<!$N%3%a%s%H$H$j$5$k(B. $B>e$N9T$r%3%a%s%H%"%&%H(B.
%     test5.data tfbToTree /rr set
     /cone.getConeInfo.rr0 rr def

     rr (DIM) getNode /cdim set
     cdim 2 get 0 get 0 get 0 get to_univNum /cdim set
     % polymake $B$N(B DIM $B$O0l$D>.$5$$$N$G(B 1 $BB-$9(B.
     cdim (1).. add /cdim set

     rr (FACETS) getNode tag 0 eq {
     % FACETS $B$r;}$C$F$$$J$$$J$i:FEY7W;;$9$k(B. 
     % POINTED, NOT__POINTED $B$bF@$i$l$k(B
       cone.debug { (Calling polymake FACETS.) message } { } ifelse
       [(FACETS) ceq2] doPolymake 1 get /rr set
       cone.debug { (Done.) message } { } ifelse
   } {  } ifelse

     rr (VERTICES) getNode tag 0 eq {
       (internal error: VERTICES is not found.) error
     } {  } ifelse

     /cone.getConeInfo.rr1 rr def
    
     rr (NOT__POINTED) getNode tag 0 eq {
       % cone $B$,(B pointed $B$N;~$O(B mmc $B$OC10L9TNs(B. genUnitMatrix $B$r;H$&(B.
       % VERTICES $B$h$j0l$D>.$5$$%5%$%:(B.
       /mmc 
         [ rr (VERTICES) getNode 2 get 0 get 0 get rest]
         genUnitMatrix 
       def
       /mmL [ ] def
       /ppt (POINTED) def
     } {
       % pointed $B$G$J$$>l9g(B,
       % cone $B$N@~7AItJ,6u4V$r7W;;(B.
       rr getLinearitySubspace /mmL set
       [mmL genUnit] mmL join /mmc set % [1,0,0,...] $B$rB-$9(B.
        mmc  asir_matrix_kernel  /mmc set % $BJd6u4V(B
        mmc removeFirstFromPolymake /mmc set   % $B$R$H$D>.$5$$%5%$%:$K(B.

       [mmL genUnit] mmL join asir_matrix_image 
        removeFirstFromPolymake /mmL set
        mmL asir_matrix_image /mmL set  % Linearity space $B$r5a$a$k(B. rm 0vector
        /ppt (NOT__POINTED) def
     } ifelse
     /arg1 [[cdim mmc mmL ppt] rr] def
  ] pop
  popVariables
  arg1
} def    


/test.put {
  /dog [(dog) [[(legs) 4] ] [1 2 3 ]] [(class) (tree)] dc def
  /man [(man) [[(legs) 2] ] [1 2 3 ]] [(class) (tree)] dc def
  /ma [(mammal) [ ] [man dog]] [(class) (tree)] dc def
  /fan [ma 1 copy] def
  ma (dog) getNode /dd set
  dd 2 get /dd2 set
  dd2 1 0 put
  ma message

  fan message
} def

/test6.data
 $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])]))$
def
% tfbToTree

/arrayToTree { [(class) (tree)] dc } def

%<
% polymake $B$h$jF@$i$l$?(B TreeObject $B$+$i(B TreeObject cone $B$r@8@.$9$k(B.
% Usages: test6.data tfbToTree newCone $B$GF0:n%F%9%H(B
%>
/test6 {
  test6.data tfbToTree /rr set
  rr newCone /rr2 set
} def

%<
% Usages: doPolymakeObj newCone
%>
/newCone {
  /arg1 set
  [/polydata /cone /facets /vertices /flipped /ineq 
   /facetsv /rr] pushVariables
  [
    /polydata arg1 def
    polydata (FACETS) getNode tag 0 eq {
      (newCone : no FACETS data.) error
    } {  } ifelse
% facets $B$OM-M}?t$N>l9g@55,2=$9$k(B.  data/test11 $B$G(B $BM-M}?t$G$k(B.
    polydata (FACETS) getNode 2 get 0 get to_univNum 
    { nnormalize_vec} map /facets set
    [[ ] ] facets join shell rest removeFirstFromPolymake /facets set
    facets length 0 eq 
    {(Internal  error. Facet data is not obtained. See OpenXM_tmp.) error} { } ifelse
% vertices $B$O(B cone $B$N>e$K$"$k$N$G@0?tG\(B OK. $B@55,$+$9$k(B.
    polydata (VERTICES) getNode 2 get 0 get to_univNum 
    { nnormalize_vec} map /vertices set
    [[ ] ] vertices join shell rest removeFirstFromPolymake /vertices set
% inequalities $B$OM-M}?t$N>l9g@55,2=$9$k(B.
    polydata (INEQUALITIES) getNode 2 get 0 get to_univNum 
    { nnormalize_vec } map /ineq set
    [[ ] ] ineq join shell rest removeFirstFromPolymake /ineq set

% nextcid, nextfid $B$r2C$($k(B.  nextcid $B$O(B nextConeId $B$NN,(B. $B$H$J$j$N(B cone $BHV9f(B.
%                           nextfid $B$O(B nextFacetId $B$NN,(B. $B$H$J$j$N(B cone $B$N(B facet
%                            $BHV9f(B. 
    [(cone) [ ] 
     [
      [(facets) [ ] facets]  arrayToTree
      [(flipped) [ ] facets length newVector null_to_zero] arrayToTree
      [(facetsv) [ ] facets vertices newCone_facetsv] arrayToTree
      [(nextcid) [ ] facets length newVector.with-1 ] arrayToTree
      [(nextfid) [ ] facets length newVector.with-1 ] arrayToTree
      [(vertices) [ ] vertices]  arrayToTree
      [(inequalities) [ ] ineq] arrayToTree
     ]
    ] arrayToTree /cone set    
    /arg1 cone def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: newCone_facetv
% facet vertices newCone_facetv
% facet $B$K$N$C$F$$$k(B vertices $B$r$9$Y$FNs5s(B.
%>
/newCone_facetv {
  /arg2 set
  /arg1 set
  [/facet /vertices] pushVariables
  [
    /facet arg1 def /vertices arg2 def
    [
      0 1 vertices length 1 sub {
         /ii set
         facet vertices ii get mul isZero
         { vertices ii get } {  } ifelse
      } for
    ]
    /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
% Usages: newCone_facetsv
% facets vertices newCone_facetv
% facets $B$K$N$C$F$$$k(B vertices $B$r$9$Y$FNs5s(B. $B%j%9%H$r:n$k(B.
%>
/newCone_facetsv {
  /arg2 set
  /arg1 set
  [/facets /vertices] pushVariables
  [
    /facets arg1 def /vertices arg2 def
    facets { vertices newCone_facetv } map
    /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
% Usages: [gb weight] newConeGB
%  gb $B$H(B weight $B$r(B tree $B7A<0$K$7$F3JG<$9$k(B.
%>
/newConeGB {
  /arg1 set
  [/gbdata  /gg /ww /rr] pushVariables
  [
    /gbdata arg1 def
% gb
    gbdata 0 get /gg set
% weight
    gbdata 1 get /ww set
% 
    [(coneGB) [ ] 
     [
      [(grobnerBasis) [ ] gg]  arrayToTree
      [(weight) [ ] [ww]] arrayToTree
      [(initial) [ ] gg { ww 2 get weightv init } map ] arrayToTree
     ]
    ] arrayToTree /rr set
    /arg1 rr def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: cone_random
%>
/cone_random.start  (2)..  def
/cone_random {
  [(tdiv_qr) 
   cone_random.start  (1103515245).. mul
   (12345).. add 

   (2147483646)..
  ] mpzext 1 get /cone_random.start set
  cone_random.start
} def

/cone_random.limit 40 def
/cone_random_vec {
  /arg1 set
  [/nn /rr] pushVariables
  [
    /nn arg1 def
    [
      0 1 nn 1 sub {
        pop
        [(tdiv_qr) cone_random  cone_random.limit] mpzext 1 get
      } for
    ] /arg1 set
  ] pop
  popVariables
  arg1 
} def

%<
% Usages: getNewRandomWeight
%%  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.
%%  h, H  $B$N=hM}$bI,MW(B. 
%% $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?
%>
/getNewRandomWeight {
  /arg1 set
  [/vv /vvd /rr] pushVariables
  [
    /vv arg1 def
    vv { (D) 2 1 roll 2 cat_n } map /vvd set
  ] pop
  popVariables
  arg1
} def

% test7 : univNum $B$N(B weight $B$,@5$7$/G'<1$5$l$k$+$N%F%9%H(B
% aux-cone.sm1

%<
% Usages: n d coneEqForSmallFan.2  (cone.type 2 $B@lMQ(B:  x,y,Dx,Dy,h)
%  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.
%  $B$O$8$a$+$i(B d $B8D$NJQ?t(B.
% 4, 2 , s,t,x,y $B$J$i(B weight $B$O(B s,t,Ds,Dt $B$N$_(B.
% u_i + v_i >= 0 ,  u_i = v_i = 0.
% homog $BJQ?t$N>r7o(B u_i+v_i >= t, i.e, -t >= 0  $B$bF~$l$k(B.
%  coneEq $B$N7k2L$H(B coneEqForSmallFan.2 $B$N7k2L$r(B join $B$7$F(B
%  getConeInfo or newCone
% note-cone.sm1  2004.8.31 $B$r8+$h(B.  w_ineq $B$"$?$j(B.
% cone.local $B$,@_Dj$5$l$F$$$k$H(B u_i <= 0 $B$b>r7o$KF~$k(B.
%>
/coneEqForSmallFan.2 {
  /arg2 set
  /arg1 set
  [/n /d /nn /dd /ii /tt] pushVariables
  [
     /n arg1 def
     /d arg2 def
     n to_int32 /n set 
     d to_int32 /d set
     /dd n d add def
     /nn n n add def

     % 0 ~ d-1, n ~ dd-1  $B$G$O(B u_i + v_i = 0
     % d ~ n-1, dd ~ nn-1 $B$G$O(B u_i=v+i = 0.
     % -t >= 0
     [
     % d ~ n-1, dd ~ nn-1 $B$G$O(B u_i=v+i = 0.
       d 1 n 1 sub {
         /ii set
      % [ 0,0, ..., 0,1,0,... ; 0] $B$r@8@.(B
         nn 1 add newVector null_to_zero  /tt set
         tt ii (1).. put
         tt    
      % [ 0,0, ..., 0,-1,0,... ; 0] $B$r@8@.(B
         nn 1 add newVector null_to_zero  /tt set
         tt ii (-1).. put
         tt
       } for
       dd 1 nn 1 sub {
         /ii set
         nn 1 add newVector null_to_zero  /tt set
         tt ii (1).. put
         tt    
         nn 1 add newVector null_to_zero  /tt set
         tt ii (-1).. put
         tt
       } for

     % 0 ~ d-1, n ~ dd-1  $B$G$O(B u_i + v_i = 0
       0 1 d 1 sub {
         /ii set
         nn 1 add newVector null_to_zero  /tt set
         tt ii (1).. put
         tt ii n add (1).. put
         tt

         nn 1 add newVector null_to_zero  /tt set
         tt ii (-1).. put
         tt ii n add (-1).. put
         tt

       } for
  
     % -t >= 0
      cone.h0 { 
      % t = 0
       nn 1 add newVector null_to_zero /tt set
       tt nn (1).. put
       tt
       nn 1 add newVector null_to_zero /tt set
       tt nn (-1).. put
       tt
      }
      {
      % -t >= 0
       nn 1 add newVector null_to_zero /tt set
       tt nn (-1).. put
       tt
      } ifelse

     % cone.local $B$,(B 1 $B$N;~(B
     % 0 ~ d-1  $B$G$O(B -u_i >= 0
      cone.local {
       0 1 d 1 sub {
         /ii set
         nn 1 add newVector null_to_zero  /tt set
         tt ii (-1).. put
         tt
       } for
      } {  } ifelse
     ] /rr set
     /arg1 rr to_univNum def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: n d coneEqForSmallFan.1  (cone.type 1 $B@lMQ(B:  x,y,Dx,Dy,h,H)
%                 cone.type 2 $B$G$O(B x,y,Dx,Dy,h 
%   coneEqForSmallFan.2 $B$N7k2L$rMQ$$$F@8@.(B.
%   H $B$N>r7o$r2C$($k(B.
%>
/coneEqForSmallFan.1 {
  /arg2 set
  /arg1 set
  [/n /d /i /j /rr /tt /tt2] pushVariables
  [
    /n arg1 def /d arg2 def
    n d coneEqForSmallFan.2 /rr set
    rr cone.appendZero /rr set
% H $BMQ$N(B 0 $B$r2C$($k(B.
% $B$H$j$"$($:(B t' = 0 $B$G$-$a$&$A(B.
    cone.h0 { } { (cone.h0 = 0 has not yet been implemented.) error } ifelse
    n 2 mul 2 add newVector null_to_zero /tt set
    tt n 2 mul 2 add 1 sub (-1).. put 
    n 2 mul 2 add newVector null_to_zero /tt2 set
    tt2 n 2 mul 2 add 1 sub (1).. put 
    rr [tt tt2] join /rr set
    /arg1 rr to_univNum def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: vv ineq toQuotientCone
% weight space $B$N(B $B%Q%i%a!<%?$D$1$N$?$a$K;H$&(B.
% cone.V $B$r5a$a$?$$(B.  vv $B$O(B doPolymakeObj (VERTICES) getNode 2 get 0 get $B$GF@$k(B.
% vertices $B$N(B non-negative combination $B$,(B cone.
% vertice cone.w_ineq isInLinearSubspace $B$J$i<h$j=|$/(B.
% $B$D$^$j(B vertice*cone.w_ineq = 0 $B$J$i<h$j=|$/(B.
%
% $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)
% cone.w_cone 1 get (VERTICES) getNode :: $B$HHf3S$;$h(B.
%  $B$3$N4X?t$r8F$s$G(B cone.W $B$r:n$k$N$OITMW$+$b(B.
%
% Example:  cf. parametrizeSmallFan
%   4 2 coneEqForSmallFan.2 /cone.w_ineq set cone.w_ineq getConeInfo /rr set 
%   rr 1 get (VERTICES) getNode 2 get 0 get removeFirstFromPolymake /vv set
%   vv cone.w_ineq toQuotientCone pmat 
%>
/toQuotientCone {
  /arg2 set /arg1 set
  [/vv /ineq /rr] pushVariables
  [
    /vv arg1 def /ineq arg2 def
    vv {
      dup
      ineq isInLinearSpace 1 eq { pop }
      {  } ifelse
    } map /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
% Usages:  n d parametrizeSmallFan
%  n : x $BJQ?t$N?t(B.
%  d : 0 $B$K$7$J$$(B weight $B$N?t(B.
% $B<!$NBg0hJQ?t$b@_Dj$5$l$k(B.
% cone.W :  weight $B$r%Q%i%a!<%?$E$1$9$k%Y%/%H%k$NAH(B.
% 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,  
%             i $B$,(B Wpos ~ $B$NHO0O$N$H$-(B V[i] $B$X$O(B Z $B$N85$r3]$1;;$7$F$h$$(B.
% cone.w_ineq :  weight space $B$NITEy<0@)Ls(B.  $B0J8e$N7W;;$G>o$KIU2C$9$k(B.
% cone.w_cone :  w_ineq $B$r(B polymake $B$G(B getConeInfo $B$7$?7k2L(B.
% Example: /cone.local 1 def ; 4 2 parametrizeSmallFan pmat
% Example: /cone.local 0 def ; 4 2 parametrizeSmallFan pmat
%>
/parametrizeSmallFan {
  /arg2 set /arg1 set
  [/n /d /vv /coneray] pushVariables
  [
    /n arg1 def /d arg2 def
    {
      cone.type 1 eq {
        n d coneEqForSmallFan.1 /cone.w_ineq set 
        exit
      } {  } ifelse
      cone.type 2 eq {
        n d coneEqForSmallFan.2 /cone.w_ineq set 
        exit
      } {  } ifelse
      (This cone.type has not yet been implemented.) error
    } loop
    cone.w_ineq getConeInfo /cone.w_cone set 
    cone.w_cone 1 get (VERTICES) getNode 2 get 0 get 
      removeFirstFromPolymake /vv set

    vv cone.w_ineq toQuotientCone  /coneray set
    coneray length /cone.Wpos set

    coneray cone.w_cone 0 get 2 get join /cone.W set
    /arg1 cone.W def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: n d coneEqForTotalFan.2  (cone.type 2 $B@lMQ(B:  x,y,Dx,Dy,h)
%  n $BJQ?t$N?t(B,
%  d 0 $B$K$7$J$$JQ?t(B.  
% u_i + v_i >= 0 , 
% homog $BJQ?t$N>r7o(B u_i+v_i >= 0, t = 0  $B$bF~$l$k(B.
%  coneEq $B$N7k2L$H(B coneEqForSmallFan.2 $B$N7k2L$r(B join $B$7$F(B
%  getConeInfo or newCone
% cone.local $B$,@_Dj$5$l$F$$$k$H(B u_i <= 0 $B$b>r7o$KF~$k(B.
%>
/coneEqForTotalFan.2 {
  /arg2 set
  /arg1 set
  [/n /nn /dd /ii /tt] pushVariables
  [
     /n arg1 def
     /d arg2 def
     n to_int32 /n set 
     d to_int32 /d set 
     /nn n n add def
     /dd n d add def

     % 0 ~ d-1, n ~ dd-1  $B$G$O(B u_i + v_i >= 0
     % d ~ n-1, dd ~ nn-1 $B$G$O(B u_i=v+i = 0.
     % t = 0
     [
     % d ~ n-1, dd ~ nn-1 $B$G$O(B u_i=v+i = 0.
       d 1 n 1 sub {
         /ii set
      % [ 0,0, ..., 0,1,0,... ; 0] $B$r@8@.(B
         nn 1 add newVector null_to_zero  /tt set
         tt ii (1).. put
         tt    
      % [ 0,0, ..., 0,-1,0,... ; 0] $B$r@8@.(B
         nn 1 add newVector null_to_zero  /tt set
         tt ii (-1).. put
         tt
       } for
       dd 1 nn 1 sub {
         /ii set
         nn 1 add newVector null_to_zero  /tt set
         tt ii (1).. put
         tt    
         nn 1 add newVector null_to_zero  /tt set
         tt ii (-1).. put
         tt
       } for

     % 0 ~ d-1, n ~ dd-1  $B$G$O(B u_i + v_i >= 0
       0 1 d 1 sub {
         /ii set
         nn 1 add newVector null_to_zero  /tt set
         tt ii (1).. put
         tt ii n add (1).. put
         tt

       } for
  
     % t = 0
      cone.h0 { 
      % t = 0
       nn 1 add newVector null_to_zero /tt set
       tt nn (1).. put
       tt
       nn 1 add newVector null_to_zero /tt set
       tt nn (-1).. put
       tt
      }
      {
         (coneForTotalFan.2. Not implemented.) error 
      } ifelse

     % cone.local $B$,(B 1 $B$N;~(B
     % 0 ~ d-1  $B$G$O(B -u_i >= 0
      cone.local {
       0 1 d 1 sub {
         /ii set
         nn 1 add newVector null_to_zero  /tt set
         tt ii (-1).. put
         tt
       } for
      } {  } ifelse
     ] /rr set
     /arg1 rr to_univNum def
  ] pop
  popVariables
  arg1
} def

%<
% Usages:  n d parametrizeTotalFan
%  n : x $BJQ?t$N?t(B.
%  d : 0 $B$K$7$J$$?t(B.
% $B<!$NBg0hJQ?t$b@_Dj$5$l$k(B.
% cone.W :  weight $B$r%Q%i%a!<%?$E$1$9$k%Y%/%H%k$NAH(B.
% 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,  
%             i $B$,(B Wpos ~ $B$NHO0O$N$H$-(B V[i] $B$X$O(B Z $B$N85$r3]$1;;$7$F$h$$(B.
% cone.w_ineq :  weight space $B$NITEy<0@)Ls(B.  $B0J8e$N7W;;$G>o$KIU2C$9$k(B.
% cone.w_ineq $B$r(B getConeInfo $B$7$?7k2L$O(B cone.w_cone
% Example: /cone.local 1 def ; 3 parametrizeSmallFan pmat
% Example: /cone.local 0 def ; 3 parametrizeSmallFan pmat
% local $B$,(B 1 $B$@$H(B u_i <= 0 $B$K$J$k(B.
%>
/parametrizeTotalFan {
  /arg2 set
  /arg1 set
  [/n /d /vv /coneray] pushVariables
  [
    /n arg1 def  /d arg2 def
    {
      cone.type 2 eq { n d coneEqForTotalFan.2 /cone.w_ineq set exit} 
      { } ifelse
      (This cone.type has not yet been implemented.) error
    } loop
    cone.w_ineq getConeInfo /cone.w_cone set 
    cone.w_cone 1 get (VERTICES) getNode 2 get 0 get 
     removeFirstFromPolymake /vv set

    vv cone.w_ineq toQuotientCone  /coneray set
    coneray length /cone.Wpos set

    coneray cone.w_cone 0 get 2 get join /cone.W set
    /arg1 cone.W def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: vlist wlist cone_wtowv
% [x y Dx Dy h] [-1 0 1 0 0] ==> [(x) -1 (Dx) 1] $B$r:n$k(B.
%>
/cone_wtowv {
  /arg2 set /arg1 set
  [/vlist /wlist /ii] pushVariables
  [
    /vlist arg1 def
    /wlist arg2 def
    wlist length vlist length eq {
    } {  (cone_wtowv: length of the argument must be the same.) error} ifelse

    wlist to_int32 /wlist set
    [
      0 1 wlist length 1 sub { 
        /ii set
        wlist ii get 0 eq { }
        { vlist ii get wlist ii get } ifelse
      } for
    ] /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
% Usages:  pruneZeroVector
%    genPo, getConeInfo $BEy$NA0$K;H$&(B.  0 $B%Y%/%H%k$O0UL#$N$J$$@)Ls$J$N$G=|$/(B.
%    $BF1$8@)Ls>r7o$b$N$>$/(B. polymake FACET $B$,@5$7$/F0$+$J$$>l9g$,$"$k$N$G(B.
%    cf. pear/OpenXM_tmp/x3y2.poly, x^3+y^2, x^2+y^3 data/test15.sm1
%>
/pruneZeroVector {
  /arg1 set
  [/mm /ii /jj /tt] pushVariables
  [
    /mm arg1 def
    mm to_univNum /mm set
    [ [ ] ] mm join shell rest uniq /mm set
    [
      0 1 mm length 1 sub {
         /ii set
         mm ii get /tt set
         {
          0 1 tt length 1 sub {
             /jj set
             tt jj get (0).. eq {  } 
             { tt exit } ifelse
          } for
          exit
         } loop
      } for
    ] /arg1 set
  ] pop
  arg1
} def

%<
% Usages: a projectIneq v ,  dim(a) = n, dim(v) = d
%  a*cone.Wt*cone.Lpt
%>
/projectIneq {
  cone.Wt mul cone.Lpt mul
} def

%<
% Usages: v liftWeight [w vw],  dim(v) = d, dim(w) = n, vw : vw $B7A<0$N(B weight
%   v*cone.Lp*cone.W   cone.vlist w cone_wtowv
%>
/liftWeight {
  /arg1 set
  [/v /w /vw] pushVariables
  [
    /v arg1 def
    v cone.Lp mul cone.W mul /w set
    [w  cone.vlist w cone_wtowv] /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
% Usage: m isZero
% dr.sm1 $B$X0\$9(B.
%>
/isZero {
  /arg1 set
  [/mm /ans /ii] pushVariables
  [
    /mm arg1 def
    /ans 1 def
    mm isArray {
      0 1 mm length 1 sub {
        /ii set
        mm ii get isZero /ans set
        ans 0 eq { exit } {  } ifelse
      } for
    } {
      {
        mm tag 1 eq {/ans mm 0 eq def exit} { } ifelse
        mm isPolynomial { /ans mm (0). eq def exit } { } ifelse
        mm isUniversalNumber { /ans mm (0).. eq def exit } { } ifelse
        /ans 0 def exit
      } loop
    } ifelse    
    /arg1 ans def
  ] pop
  popVariables
  arg1
} def
[(isZero)
[(m isZero bool)]] putUsages

%<
% Usage: m isNonNegative
% dr.sm1 $B$X0\$9(B.
%>
/isNonNegative {
  /arg1 set
  [/mm /ans /ii] pushVariables
  [
    /mm arg1 def
    /ans 1 def
    mm isArray {
      0 1 mm length 1 sub {
        /ii set
        mm ii get isNonNegative /ans set
        ans 0 eq { exit } {  } ifelse
      } for
    } {
      {
        mm tag 1 eq {/ans mm 0 gt mm 0 eq or def exit} { } ifelse
        mm isUniversalNumber { /ans mm (0).. gt mm (0).. eq or def exit } 
        { } ifelse
        mm isRational { mm (numerator) dc mm (denominator) dc mul /mm set
          /ans mm (0).. gt mm (0).. eq or def exit } {  } ifelse
        /ans 0 def exit
      } loop
    } ifelse    
    /arg1 ans def
  ] pop
  popVariables
  arg1
} def
[(isNonNegative)
[(m isNonNegative bool)
 (In case of matrix, m[i,j] >= 0 must hold for all i,j.)
]] putUsages

% Global variable:  cone.weightBorder
% /cone.weightBorder null def  $BITMW$G$"$m$&(B.  getStartingCone $B$G@_Dj$5$l$k(B.

%<
% Usages: cone i isOnWeigthBorder
% cone $B$N(B i $BHVL\$N(B facet $B$,(B weight $B6u4V$N6-3&$K$"$k$+(B?
% $BBg0hJQ?t(B  cone.weightBorder $B$,@_Dj$5$l$F$k$3$H(B.
% $B$3$NJQ?t$O(B cone $B$N(B facet $B%Y%/%H%k$N%j%9%H(B.
% $B$3$NJQ?t$O(B setWeightBorder $B$G@_Dj(B
% cone.weightBorder[0] or cone.weightBorder[1] or ...
% /ccone cone.startingCone def  ccone 0 isOnWeightBorder
%                               ccone 1 isOnWeightBorder
%> 
/isOnWeightBorder {
  /arg2 set /arg1 set
  [/cone /facet_i /i /j /vv /co /ans] pushVariables
  [
    /cone arg1 def /facet_i arg2 def
    facet_i to_int32 /facet_i set
    /ans 0 def
    cone (facetsv) getNode 2 get facet_i get /vv set % Facet $B$r(B vertex $BI=8=(B.
    {
      0 1 cone.weightBorder length 1 sub {
         /i set
         cone.weightBorder i get /co set % co $B$K@)Ls>r7o(B
         vv cone.Lp mul  % vv $B$r(B weight space $B$X(B lift.
         co mul isZero 
         { /ans 1 def exit }  {   } ifelse
      } for
      exit
    } loop 
    /arg1 ans def
  ] pop 
  popVariables
  arg1
} def

%<
% Usages: cone i markFlipped
% 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.
% cone $B$O(B class-tree.  Constructor $B$O(B newCone
%> 
/markFlipped {
  /arg2 set /arg1 set
  [/cone /facet_i /vv] pushVariables
  [
    /cone arg1 def /facet_i arg2 def
    facet_i to_int32 /facet_i set
    cone (flipped) getNode 2 get /vv set
    vv facet_i (1).. put
  ] pop 
  popVariables
} def

%<
% Usages: cone i [cid fid] markNext
% 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.
%   cone $B$N(B nextcid[i] = cid; nextfid[i] = fid $B$H$J$k(B.
% cone $B<+BN$,JQ99$5$l$k(B.
% cone $B$O(B class-tree.
%> 
/markNext {
  /arg3 set /arg2 set /arg1 set
  [/cone /facet_i /vv /nextid] pushVariables
  [
    /cone arg1 def /facet_i arg2 def /nextid arg3 def
    facet_i to_int32 /facet_i set
    cone (nextcid) getNode 2 get /vv set
    vv facet_i , nextid 0 get to_univNum , put

    cone (nextfid) getNode 2 get /vv set
    vv facet_i , nextid 1 get to_univNum , put
  ] pop 
  popVariables
} def



%<
% Usages: cone getNextFacet i
% flipped $B$N(B mark $B$N$J$$(B facet $B$N(B index facet_i $B$rLa$9(B.
% $B$=$l$,$J$$$H$-$O(B null
%> 
/getNextFacet {
  /arg1 set
  [/cone /facet_i /vv /ii] pushVariables
  [
    /cone arg1 def 
    /facet_i null def
    cone (flipped) getNode 2 get /vv set
    0 1 vv length 1 sub {
       /ii set
       vv ii get to_int32 0 eq { /facet_i ii def exit } 
       {  } ifelse
    } for
    /arg1 facet_i def
  ] pop 
  popVariables
  arg1
} def

%<
% Usages: cone i epsilon flipWeight
% cone $B$N(B i $BHVL\$N(B facet $B$K$+$s$7$F(B flip $B$9$k(B.
% $B?7$7$$(B weight $B$r5a$a$k(B.  cf. liftWeight
%> 
/flipWeight {
  /arg3 set /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
    /ep arg3 def

    ep to_univNum (1).. div /ep set 

% note: 2004.9.2
    cone (facetsv) getNode 2 get facet_i get /v set
    cone (facets)  getNode 2 get facet_i get /f set
    /vp v 0 get def
    1 1 v length 1 sub {
      /ii set
      vp v ii get  add /vp set
    } for
    vp ep f mul sub /vp set
    vp nnormalize_vec /vp set
    /arg1 vp def
  ] pop 
  popVariables
  arg1 
} def

%<
% Usages: cone1 cone2 isSameCone bool
% cone1 cone2 $B$,Ey$7$$$+(B? facet $B$GHf$Y$k(B.
% cone1, cone2 $B$O(B pointed cone $B$G$J$$$H$$$1$J$$(B.
%> 
/isSameCone {
  /arg2 set /arg1 set
  [/cone1 /cone2 /facets1 /facets2 /ans] pushVariables
  [
    /cone1 arg1 def 
    /cone2 arg2 def
    /facets1  cone1 (facets) getNode 2 get def
    /facets2  cone2 (facets) getNode 2 get def
    facets1 length facets2 length eq {
      facets1 facets2 sub isZero /ans set
    } {  
      /ans 0 def
    } ifelse
    /arg1 ans def 
  ] pop 
  popVariables
  arg1
} def

%<
% Usages: cone1 cone2 getCommonFacet list
% cone1 $B$NCf$G(B cone2 $B$K4^$^$l$k(B facet $B$N%j%9%H(B
% cone2 $B$NCf$G(B cone1 $B$K4^$^$l$k(B facet $B$N%j%9%H$r$b$I$9(B.
%  [1 [i] [j]] $B$"$k$H$-(B.  [0 [ ] [ ]] $B$J$$$H$-(B.
% cone1 $B$N(B facetsv[i] $B$,(B cone2 $B$K4^$^$l$k$+D4$Y$k(B.
% cone2 $B$N(B facetsv[i] $B$,(B cone1 $B$K4^$^$l$k$+D4$Y$k(B.
% cone1, cone2 $B$O(B pointed cone $B$G$J$$$H$$$1$J$$(B.
%> 
/getCommonFacet {
  /arg2 set /arg1 set
  [/cone1 /cone2 /facets /ineq /ans1 /ans2 /i /tt] pushVariables
  [
    /cone1 arg1 def 
    /cone2 arg2 def

    /facets  cone1 (facetsv) getNode 2 get def
    /ineq cone2 (inequalities) getNode 2 get def
    /ans1 [
      0 1 facets length 1 sub {
        /i set
        facets i get /tt set % facetsv[i] $B$r(B tt $B$X(B.
        ineq tt transpose mul isNonNegative {
          i
        } {  } ifelse
      } for
    ] def

    /facets  cone2 (facetsv) getNode 2 get def
    /ineq cone1 (inequalities) getNode 2 get def
    /ans2 [
      0 1 facets length 1 sub {
        /i set
        facets i get /tt set % facetsv[i] $B$r(B tt $B$X(B.
        ineq tt transpose mul isNonNegative {
          i
        } {  } ifelse
      } for
    ] def
    ans1 length 1 gt ans2 length 1 gt or {
      (getCommonFacet found more than 1 common facets.) error
    } {  } ifelse
% $B6&DL(B facet $B$,$"$l$P(B 1, $B$J$1$l$P(B 0.
    ans1 length 1 eq ans2 length 1 eq and {
      /tt 1 def
    } {
      /tt 0 def
    } ifelse
    /arg1 [tt ans1 ans2] def
  ] pop 
  popVariables
  arg1
} def

%
% -------------------------------------------------
% test8 $B$O(B aux-cone.sm1 $B$X0\F0(B.
% $B0J2<$$$h$$$h0lHL$N%W%m%0%i%`$N:n@.3+;O(B.
% -------------------------------------------------
%

%<
% Usages: setWeightBorder
%  cone.weightBorder (weight cone $B$N(B facet $B%Y%/%H%k$N=89g(B) $B$r@_Dj$9$k(B.
%  $B$"$HI{;:J*$H$7$F(B  cone.w_cone_projectedWt (doPolymakeObj)
%                    cone.w_ineq_projectedWt
%  cone.m $B<!85$N%Y%/%H%k(B.
%  cone.W, cone.Wt, cone.w_ineq $B$,$9$G$K7W;;$:$_$G$J$$$H$$$1$J$$(B.
%> 
/setWeightBorder {
  [
    (Entering setWeightBorder ) message
    cone.w_ineq cone.Wt mul pruneZeroVector /cone.w_ineq_projectedWt set
    {
      cone.w_ineq_projectedWt length 0 eq { 
% weight $B$N6u4V$K(B border $B$,$J$$>l9g(B.
        /cone.weightBorder [ ] def
        exit
      } {  } ifelse
% weight $B$N6u4V$K(B border $B$,$"$k>l9g(B.
      cone.w_ineq_projectedWt getConeInfo /cone.w_cone_projectedWt set
      cone.w_cone_projectedWt 0 get 0 get to_int32 cone.m to_int32 eq {
      } {
        (setWeightBorder : internal error.) message
      } ifelse
      cone.w_cone_projectedWt 1 get (FACETS) getNode 2 get 0 get
      removeFirstFromPolymake /cone.weightBorder set
      exit
    } loop
    (cone.weightBorder=) message
    cone.weightBorder pmat
  ] pop
} def

%
% -------------------------------------------------
% $B%W%m%0%i%`$NN.$l(B.
% Global: cone.fan   cone $B$rG[Ns$H$7$F3JG<$9$k(B.
%
% ncone (next cone) $B$,?75,$KF@$i$l$?(B cone $B$G$"$k$H$9$k(B.
% $B$3$N$H$-<!$NA`:n$r$9$k(B.
%  0. ncone $B$,(B cone.fan $B$K$9$G$K$J$$$+D4$Y$k(B. $B$"$l$P(B, internal error.
%  1. ncone markBorder ; ncone $B$NCf$N(B border $B>e$N(B facet $B$r(B mark
%  2. cone.fan $B$NCf$N(B cone $B$H6&DL(B facet $B$,$J$$$+D4$Y(B (getCommonFacet),
%     $B$"$l$P$=$l$i$r(B mark $B$9$k(B.
%     global: cone.incidence $B$K(B $B6&DL(Bfacet $B$r;}$DAH$_$N>pJs$r2C$($k(B.
%  3. ncone $B$r(B cone.fan $B$N:G8e$K2C$($k(B.
%  $B0J>e$NA`:n$r$^$H$a$?$b$N$,(B  ncone updateFan
%
%  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.
%  $B$J$1$l$P(B null $B$rLa$9(B.  null $B$,La$l$P%W%m%0%i%`=*N;(B.
%
%  getStargingCone $B$O7W;;$r=PH/$9$Y$-?75,$N(B cone $B$r7W;;$9$k(B. $BBg0hJQ?t(B cone.Lt, cone.W
%  $B$J$I$b$3$NCf$G@_Dj$9$k(B.
%  $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
%  $B$H$7$FF~NO$7$F$*$/(B.
%
%  reduced gb $B$O(B $B4X?t(B input weight cone.gb reduced_G $B$G7W;;$9$k(B.
%
%
%  [ccone i] getNextCone ncone : flip $B$K$h$j<!$N(B cone $B$rF@$k(B.
%
%  1. clearGlobals ; $BF~NOBg0hJQ?t$N@_Dj(B.
%  2. getStartingCone /ncone set  
%  3. {  ncone updateFan
%  4.    getNextFlip /cone.nextflip set
%  6.    cone.nextflip isNull { exit } {  } ifelse
%  7.    cone.nextflip getNextCone /ncone set
%  8. } loop
%  
%
% -------------------------------------------------
%

%<
% Usages: input weight cone.gb_Dh reduced_G
%  gb in h[1,1](D)
%>
/cone.gb_Dh {
  /arg2 set /arg1 set
  [/ff /ww /gg] 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
    /cone.gb_Dh.g gg def
    /arg1 gg def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: cone.boundp
%
/cone.boundp {
   dup boundp 2 1 roll tag 0 eq not and
} def

%<
% Usages: clearGlobals
% cf. cone.boundp
% polymake $B$r:FEY8F$V$?$a$K(B global $BJQ?t$r%/%j%"$9$k(B.
% $B$^$@ESCf(B.
%>
/clearGlobals {
  /cone.W null def
  /cone.Wt null def

  /cone.cinit null def
  /cone.weightBorder null def

} def

%<
% Usages: getStartingCone ncone
% getStargingCone $B$O7W;;$r=PH/$9$Y$-?75,$N(B cone $B$r7W;;$9$k(B.
% $B@_Dj$9$Y$-Bg0hJQ?t$O0J2<$r8+$h(B.
%>

/getStartingCone.test {
%------------------Globals----------------------------------------
% ---------------  $BF~NO%G!<%?MQBg0hJQ?t$N@_Dj(B --------------------------
%
% cone.input : $BF~NOB?9`<07O(B
/cone.input 
  [(t1-x-y) (h*t2-x^2-y^2) (2*x*Dt2+h*Dt1+h*Dx) (2*y*Dt2+h*Dt1+h*Dy)]
def

% cone.vlist : $BA4JQ?t$N%j%9%H(B
/cone.vlist [(t1) (t2) (x) (y) (Dt1) (Dt2) (Dx) (Dy) (h)] def

% cone.vv : define_ring $B7A<0$NJQ?t%j%9%H(B.
% t1,t2, x,y   : t-space $B$N(B Grobner fan (local) $B$r5a$a$k(B.
/cone.vv (t1,t2,x,y) def

% cone.parametrizeWeightSpace : weight $B6u4V$r(B parametrize $B$9$k4X?t(B. 
%   $BBg0hJQ?t(B cone.W , cone.Wpos $B$b$-$^$k(B.
/cone.parametrizeWeightSpace {
  4 2 parametrizeSmallFan
} def

% cone.w_start : weight$B6u4V$K$*$1$k(B weight $B$N=i4|CM(B. 
% $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. 
/cone.w_start
  [ 1 4 ]
def   

% cone.gb : gb $B$r7W;;$9$k4X?t(B.
/cone.gb {
  cone.gb_Dh
} def

%
% -----------------  $B$*$o$j(B ---------------------------
%
} def  % end of getStartingCone.test

/getStartingCone {
 [/wv_start /w_start /reduced_G] pushVariables
 [
% cone.n $B$O<+F0E*$K$-$a$i$l$k(B.
%  cone.n $B$O(B GB $B$r7W;;$9$k6u4V$N<!85(B.
  /cone.n cone.vlist length def
%[1]  cone.W, cone.Wpos $B$r5a$a$k(B.   cone.m $B$O(B cone.W $B$h$j<+F0E*$K$-$^$k(B.
%  cone.m $B$O(B weight $B6u4V$N<+M3EY(B. cone.W $B$G<M1F$5$l$k@h$N<!85(B.
  /cone.W cone.boundp {
    (Skip cone.parametrizeWeightSpace. cf. clearGlobals) message
  } {
    cone.parametrizeWeightSpace
  } ifelse
  (parametrizing weight space: cone.W = ) messagen cone.W message
  /cone.Wt cone.W transpose def
  /cone.m cone.W length def
% WeightBorder $B$N>r7oH=Dj(B facet $B$r@_Dj(B.
  /cone.weightBorder cone.boundp {
    (Skip setWeightBorder cf. clearGlobals) message
  } {
    setWeightBorder
  } ifelse

%[2] weight vector wv_start $B$r@8@.$9$k(B.
% wv_start $B$r@_Dj(B.
  cone.w_start tag 0 eq {
% cone.w_start $B$,(B null $B$J$i(B random $B$K(B weight $B$r@_Dj(B.
    /cone.w_start cone.m cone_random_vec def
  } {  
    cone.w_start length cone.m to_int32 eq {
    } {
      (Error: cone.w_start has wrong length.) error
      /cone.w_start cone.m cone_random_vec def
    } ifelse
  } ifelse 
  /w_start cone.w_start cone.W mul def 

  {
     cone.vlist w_start cone_wtowv /wv_start set
     (Trying a starting weight vector : ) messagen
     wv_start pmat
%[3] reduced GB $B$N7W;;(B.
     cone.input wv_start cone.gb /reduced_G set
     (Reduced GB is obtained: ) message
     %reduced_G pmat 
     /cone.cgb reduced_G def
     [cone.w_start w_start wv_start] /cone.cgb_weight set

%[4] $B<M1F$7$F$+$i(B polytope $B$N%G!<%?$r7W;;(B.  
     wv_start reduced_G coneEq /cone.g_ineq set
     cone.g_ineq cone.w_ineq join  /cone.gw_ineq set
     cone.gw_ineq  cone.Wt mul /cone.gw_ineq_projectedWt set % $B<M1F(B
     /cone.cinit cone.boundp {
       (Skipping cone.gw_ineq_projectedWt getConeInfo. cf. clearGlobals) message
     } {
      cone.gw_ineq_projectedWt getConeInfo /cone.cinit set
     } ifelse

     (cone.cinit is --- the first number is the dim of cone.) messagen
     cone.cinit 0 get pmat
% Maximal dimensional cone $B$+$I$&$+$N8!::(B. $B8!::$K%Q%9$9$l$P(B loop $B$r(B exit
% $B%Q%9$7$J$$>l9g(B  w_start $B$r(B cone_random_vec $B$rMQ$$$FJQ99$9$k(B.
     cone.cinit 0 get 0 get to_int32 cone.m eq { exit }
     {
       (Failed to get the max dim cone. Updating the weight ...) messagen
       cone.m cone_random_vec /cone.w_start set
       /w_start  cone.w_start cone.W mul def
% cone.cinit $B$r:FEY7W;;$9$k$?$a$K(B clear $B$9$k(B.
       /cone.cinit null def
     } ifelse
  } loop

  (cone.m = ) messagen cone.m message
  (Suceeded to get the maximal dimensional startingCone.) message

% Linearity subspace $B$N(B orth complement $B$X$N<M1F9TNs(B. 
% $BBg0hJQ?t(B cone.Lp, cone.Lpt $B$r@_Dj(B
  cone.cinit 0 get 1 get /cone.Lp set
  cone.Lp transpose /cone.Lpt set
% Linearity subspace $B$N9TNs$r@_Dj(B. 
% $BBg0hJQ?t(B cone.L $B$r@_Dj(B
  cone.cinit 0 get 2 get /cone.L set
% 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.
% $BBg0hJQ?t(B cone.d $B$N@_Dj(B.
  /cone.d cone.Lp length def

  cone.m cone.d  eq {
    (There is no linearity space) message
  } {  
    (Dim of the linearity space is ) messagen cone.m cone.d sub message
    (cone.Lp = ) messagen cone.Lp pmat 
  } ifelse

%[5] cone.g_ineq * cone.Wt * cone.Lpt 
%    cone.w_ineq * cone.Wt * cone.Lpt
%   $B$G@)Ls$r(B d $B<!85%Y%/%H%k$KJQ49(B.
% W (R^m) $B6u4V$NITEy<0@)Ls$r(B L' (R^d) $B6u4V$X<M1F(B
% cone.gw_ineq_projectedWtLpt 
%  = cone.g_ineq*cone.Wt*cone.Lpt \/ cone.w_ineq*coneWt*cone.Lpt

  /cone.gw_ineq_projectedWtLpt 
     cone.gw_ineq_projectedWt cone.Lpt mul 
  def

  cone.m cone.d eq  {
    /cone.cinit.d cone.cinit def
  } {
% cone.m > cone.d $B$J$i$P(B, $B:FEY(B cone $B$N7W;;$,I,MW(B.
% R^d $B$N(B cone $B$O(B cone.cinit.d $B$XF~$l$k(B.
    cone.gw_ineq_projectedWtLpt getConeInfo /cone.cinit.d set
  } ifelse

  cone.cinit.d 1 get newCone /cone.startingCone set

  (cone.startingCone is ) message
  cone.startingCone message
 ] pop
 popVariables
 cone.startingCone
} def

%
%  data/test9.sm1 $B$N(B test9   1-simplex X 2-simplex
% 
%  data/test10.sm1   1-simplex X 3-simplex
%  data/test11.sm1   SST, p.59
%
%  $B$$$h$$$h(B, cone enumeration $B$N%W%m%0%i%`=q$-3+;O(B
% 

%<
% Usages: cone markBorder
%   cone->facets[i] $B$,(B weight space $B$N(B border $B$K$"$k$H$-(B    
%   cone->flipped[i] = 2 $B$H$9$k(B.
%   $B$3$l$r(B cone $B$N$9$Y$F$N(B facet $B$KBP$7$F7W;;(B.
%>
/markBorder {
  /arg1 set
  [/cone /facets_t /flipped_t /kk /nextcid_t /nextfid_t] pushVariables
  [
    /cone arg1 def
    cone (facets) getNode 2 get /facets_t set
    cone (flipped) getNode 2 get /flipped_t set
    cone (nextcid) getNode 2 get /nextcid_t set
    cone (nextfid) getNode 2 get /nextfid_t set
    0 1 flipped_t length 1 sub {
      /kk set
      flipped_t kk get (0).. eq {
         cone kk isOnWeightBorder {
% Border $B$N>e$K$"$k$N$G(B flip $B:Q$N%^!<%/$r$D$1$k(B.
           flipped_t kk (2).. put
% $B$H$J$j$N(B cone $B$N(B id (nextcid, nextfid) $B$O(B -2 $B$H$9$k(B.
           nextcid_t kk (-2).. put
           nextfid_t kk (-2).. put
         } {  } ifelse
      } {  } ifelse
    } for
  ] pop
  popVariables
} def

%<
% Usages: ncone updateFan
% $B%0%m!<%P%kJQ?t(B cone.fan $B$r99?7$9$k(B.
%>
%
% updateFan $B$N(B debug $B$O(B data/test8 $B$G$H$j$"$($:$d$k(B.
%  test8 /ncone set $B$r<B9T$7$F$+$i(B  ncone updateFan

% global: cone.fan
/cone.fan [  ] def
% global: cone.incidence
/cone.incidence [ ] def
% global: cone.gblist   gb's standing for each cones in cone.fan.
/cone.gblist [ ] def

/updateFan {
  /arg1 set
  [/ncone /kk /cfacet /ii /jj /tcone /flipped_t] pushVariables
  [
    /ncone arg1 def
    /cone.fan.n  cone.fan length def
% -1.  cone.cgb ($BD>A0$K7W;;$5$l$?(B gb) $B$H(B cone.cgb_weight ($BD>A0$N7W;;$N(B weight)
%    $B$r(B cone.gblist $B$X3JG<$9$k(B.
    cone.gblist [ [cone.cgb cone.cgb_weight] newConeGB ] join /cone.gblist set
% 0. ncone $B$,(B cone.fan $B$K$9$G$K$"$l$P%(%i!<(B
    0 1 cone.fan.n 1 sub {
      /kk set
      ncone cone.fan kk get isSameCone { 
         (Internal error updateFan: ncone is already in cone.fan) error
      } {  } ifelse
    } for

% 1. ncone $B$NCf$N(B border $B>e$N(B facet $B$r$9$Y$F(B mark.
    ncone markBorder

% 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
    0 1 cone.fan.n 1 sub {
      /kk set
      ncone cone.fan kk get getCommonFacet  /cfacet set
      cfacet 0 get 
      {
% $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. 
         /ii cfacet 1 get  0 get def
         /jj cfacet 2 get  0 get def
         cone.incidence [ [[cone.fan.n ii] [kk jj]] ] join /cone.incidence set
% flipped $B$r(B mark $B$9$k(B.
         ncone ii markFlipped
         cone.fan kk get /tcone set
         tcone jj markFlipped
% nextcid, nextfid $B$r@_Dj$9$k(B.
         ncone ii [kk jj] markNext
         tcone jj [cone.fan.n ii] markNext
      } {  } ifelse
    } for
% 3. ncone $B$r2C$($k(B.
    cone.fan [ncone] join /cone.fan set
  ] pop
  popVariables
} def

%<
% usages: getNextFlip [cone, k]
% 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.
% $B$b$&$J$$$H$-$K$O(B null $B$rLa$9(B.
%>
/getNextFlip {
  [/tcone /ans /ii ] pushVariables
  [ 
    /ans null def
    0 1 cone.fan length 1 sub {
      /ii set
      cone.fan  ii get /tcone set
      tcone getNextFacet /ans set
      ans tag 0 eq { } { exit } ifelse
    } for
    ans tag 0 eq { /arg1 null def }
    { /arg1 [tcone ans] def } ifelse
  ] pop 
  popVariables
  arg1 
} def

% global variable : cone.epsilon , cone.epsilon.limit
%   flip $B$N;~$N(B epsilon
/cone.epsilon (1).. (10).. div def
/cone.epsilon.limit (1).. (100).. div def

%<
%  Usages: result_getNextFlip getNextCone ncone
%  flip $B$7$F?7$7$$(B ncone $B$rF@$k(B.
%>
/getNextCone {
 /arg1 set
 [/ncone /ccone /kk /w /next_weight_w_wv] pushVariables
 [ 
  /ccone arg1 def
  /ncone null def
  /kk ccone 1 get def
  ccone 0 get /ccone set
  {
   ccone tag 0 eq { exit } {  } ifelse

% ccone $B$N(B kk $BHVL\$N(B facet $B$K$D$$$F(B flip $B$9$k(B.
   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.input next_weight_w_wv 1 get cone.gb /cone.cgb set
   [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
% $B<!85$rD4$Y$k(B.  $B$@$a$J$i(B 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 $B$,$J$1$l$P(B $B$d$O$j(B epsilon $B$r>.$5$/(B.
   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
   (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

%<
% Usages: set globals and getGrobnerFan
%  cf. clearGlobals
% getStartingCone $B$9$k$H(B weightSpace $B$H$+$N7W;;$,$G$-$k(B. isOnWeightBorder $B$,(B
%  $B7h$a$i$l$k(B.
%>
% $B$H$j$"$($:(B (data/test8.sm1) run $B$7$F$+$i(B getGrobnerFan
/getGrobnerFan {
  getStartingCone /cone.ncone set
  {
    cone.ncone updateFan
    (  ) message
    (----------------------------------------------------------) message
    (getGrobnerFan #cone.fan=) messagen cone.fan length message
    cone.ncone /cone.ccone set
    getNextFlip /cone.nextflip set
    cone.nextflip tag 0 eq { exit } { } ifelse
    cone.nextflip getNextCone /cone.ncone set
  } loop
  (Construction  is completed. See cone.fan, cone.incidence and cone.gblist.) 
  message
} def

%<
% Usages: vlist generateD1_1
%  -1,1  weight $B$r@8@.$9$k(B.
%  vlist $B$O(B (t,x,y) $B$+(B [(t) (x) (y)]
%  
%>
/generateD1_1 {
  /arg1 set
  [/vlist /rr /rr /ii /vv] pushVariables 
  [
    /vlist arg1 def
    vlist isString {
      [vlist to_records pop] /vlist set
    } {  } ifelse
    [
      0 1 vlist length 1 sub {
        /ii set
        vlist ii get /vv set
        vv -1
        [@@@.Dsymbol vv] cat 1
      } for
    ] /rr set
    /arg1 rr def    
  ] pop
  popVariables
  arg1
} def

/listNodes {
  /arg1 set
  [/in-listNodes /ob /rr /rr /ii] pushVariables
  [
    /ob arg1 def
    /rr [ ] def
    {
      ob isClass {
        ob (array) dc /ob set
      } { exit } ifelse
      rr [ob 0 get] join /rr set
      ob 2 get /ob set
      0 1 ob length 1 sub {
         /ii set
         rr ob ii get listNodes join /rr set
      } for 
      exit
    } loop 
    /arg1 rr def
  ] pop
  popVariables
  arg1
} def
[(listNodes)
[(ob listNodes)
 (cf. getNode)
 (Example:)
 (  /dog [(dog) [[(legs) 4] ] [ ]] [(class) (tree)] dc def)
 (  /man [(man) [[(legs) 2] ] [ ]] [(class) (tree)] dc def)
 (  /ma [(mammal) [ ] [man dog]] [(class) (tree)] dc def)
 (  ma listNodes )
]] putUsages 

%<
% Usages:  obj printTree
%>
/printTree {
  /arg1 set
  [/ob /rr /rr /ii /keys /tt] pushVariables
  [
    /ob arg1 def
    /rr [ ] def
    /keys ob listNodes def
    keys 0 get /tt set
    keys rest /keys set
    keys { ob 2 1 roll getNode } map /rr set
    (begin ) messagen  tt messagen 
    ( ---------------------------------------) message
    0 1 rr length 1 sub {
       /ii set
       keys ii get messagen (=) message
       rr ii get 2 get pmat
    } for 
    (--------------------------------------- end ) messagen
    tt message
    /arg1 rr def
  ] pop
  popVariables
  arg1
} def

%<
% Usages $B$O(B (inputForm) usages $B$r$_$h(B.
%>
/inputForm {
  /arg1 set
  [/ob /rr /i ] pushVariables
  [
    /ob  arg1 def
    /rr [ ] def
    {
     ob isArray {
       rr [ ([) ] join /rr set
       0 1 ob length 1 sub {
         /i set
         i ob length 1 sub lt {
           rr [ob i get inputForm $ , $] join /rr set
         } {
           rr [ob i get inputForm] join /rr set
         } ifelse
       } for
       rr [ (]) ] join cat /rr set 
       exit 
     } { } ifelse
     ob isClass {
       ob etag 263 eq { % tree
         /rr ob inputForm.tree def exit
       } { /rr [( $ this etag is not implemented $ )] cat def exit  } ifelse
     } {  } ifelse
     ob isUniversalNumber {
       [$($ ob toString $)..$] cat /rr set 
       exit
     } {  } ifelse
     ob isPolynomial {
       [$($ ob toString $).$] cat /rr set 
       exit
     } {  } ifelse
     ob isRational {
       [$ $ ob (numerator) dc inputForm $ $
            ob (denominator) dc inputForm $ div $ ] cat /rr set
       exit
     } {  } ifelse
     ob isString {
       [$($ ob $)$ ] cat /rr set
       exit
     } {  } ifelse
     ob toString /rr set
     exit
    } loop
    rr /arg1 set
  ] pop
  popVariables
  arg1
} def
[(inputForm)
 [(obj inputForm str)
]] putUsages
% should be moved to dr.sm1 

/inputForm.tree {
  /arg1 set
  [/ob /key /rr /rr /ii] pushVariables
  [
    /ob arg1 def
    /rr [ ] def
    {
      ob (array) dc /ob set
      /rr [ $[$ ob 0 get inputForm $ , $
            ob 1 get inputForm $ , $
          ] def
      rr  [ob 2 get inputForm ] join /rr set
      rr [$ ] $] join /rr set
      rr [ $ [(class) (tree)] dc $ ] join /rr set
      rr cat /rr set
      exit
    } loop 
    /arg1 rr def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: str inputForm.value str
%>
/inputForm.value {
  /arg1 set
  [/key /val /valstr /rr] pushVariables
  [
    arg1 /key set
    key isString { } {(inputForm.value: argument must be a string) error } ifelse
    key boundp {
     [(parse) key] extension pop
     /val set
     val inputForm /valstr set
     [( ) valstr ( /) key ( set )] cat /rr set
    } {
     /valstr [] cat /rr set
    } ifelse
    rr /arg1 set
  ] pop
  popVariables
  arg1
} def

% global: cone.withGblist
/cone.withGblist 0 def
%<
% Usages:  saveGrobnerFan  str
%  GrobnerFan $B$N%G!<%?$r(B inputForm $B$KJQ99$7$FJ8;zNs$KJQ$($k(B.
%  $B$3$N%G!<%?$r(B parse $B$9$k$H(B GrobnerFan $B$rF@$k$3$H$,2DG=(B.
%  BUG: $BB?9`<0$NB0$9$k4D$N%G!<%?$NJ]B8$O$^$@$7$F$J$$(B.
%>
/saveGrobnerFan {
  [/rr] pushVariables
  [
    (cone.withGblist=) messagen cone.withGblist message
    [ 
% $B%f!<%6$N@_Dj$9$k%Q%i%a!<%?(B. cone.gb, cone.parametrizeWeightSpace $BEy$N4X?t$b$"$j(B.
      (cone.comment)
      (cone.type)  (cone.local) (cone.h0)
      (cone.vlist) (cone.vv)
      (cone.input)

% $B%W%m%0%i%`Cf$GMxMQ$9$k(B, $BBg;v$JBg0hJQ?t(B.  weight vector $B$N<M1F9TNs$,=EMW(B.
      (cone.n) (cone.m) (cone.d) 
      (cone.W) (cone.Wpos) (cone.Wt)
      (cone.L) (cone.Lp) (cone.Lpt)
      (cone.weightBorder)
      (cone.w_ineq)
      (cone.w_ineq_projectedWt)
      (cone.epsilon)

% $B7k2L$NMWLs(B.
      (cone.fan)
      cone.withGblist { (cone.gblist) } {  } ifelse
      (cone.incidence)

    ] { inputForm.value  nl } map /rr set
    rr cat /rr set
% ring $B$r(B save $B$7$F$J$$$N$GEv:B$NBP=h(B.
    [ ([) cone.vv inputForm ( ring_of_differential_operators 0 ] define_ring )
      nl nl rr] cat /arg1 set
  ] pop
  popVariables
  arg1
} def

/printGrobnerFan.1 {
  /arg1 set
  [/key /rr] pushVariables
  [
    /key arg1 def
    key boundp {
      [(parse) key] extension pop /rr set
      rr isArray {
        key messagen ( = ) message  rr pmat
      } {
        key messagen ( = ) messagen rr message
      } ifelse      
    }{  
      key  messagen ( = ) message
    } ifelse
  ] pop
  popVariables
} def

/printGrobnerFan {
  [/i] pushVariables
  [
  (==========  Grobner Fan ====================) message
   [
      (cone.comment)
      (cone.vlist) (cone.vv)
      (cone.input)
      (cone.type)  (cone.local) (cone.h0)
      (cone.n) (cone.m) (cone.d) 
      (cone.W) (cone.Wpos) (cone.Wt)
      (cone.L) (cone.Lp) (cone.Lpt)
      (cone.weightBorder)
      (cone.incidence)
   ] { printGrobnerFan.1 } map
   (   ) message
   0 1 cone.fan length 1 sub {
     /ii set
     ii messagen ( : ) messagen
     cone.fan ii get printTree
   } for
   cone.withGblist {
    0 1 cone.gblist length 1 sub {
      /ii set
      ii messagen ( : ) messagen
      cone.gblist ii get printTree
    } for
  } {  } ifelse


  (=========================================) message
  (cone.withGblist = ) messagen cone.withGblist message
  (  ) message
  ] pop
  popVariables
} def

%<
% Usages:  m uniq
% Remove duplicated lines.
%> 
/uniq  {
  /arg1 set
  [/mm /prev /i /rr] pushVariables
  [
    /mm arg1 def
   {
     mm length 0 eq { [ ] /rr set exit } {  } ifelse
     /prev mm 0 get def
     [
       prev
       1 1 mm length 1 sub {
         /i set
         mm i get prev sub isZero { } 
         { /prev mm i get def prev } ifelse
       } for
      ] /rr set
      exit
    } loop
    rr /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
% Usages: [vlist vw_vector] getGrRing [vlist vGlobal sublist]
%      example:  [(x,y,z) [(x) -1 (Dx) 1 (y) 1 (Dy) 2]] getGrRing 
%                [(x,y,z,y') [(x)] [[(Dy) (y')]]]
%  h[0,1](D_0) $B@lMQ$N(B getGrRing.
%     u_i + v_i > 0 $B$J$i(B  Dx_i ==> x_i' ($B2D49$JJQ?t(B). sublist $B$X(B.
%     u_i < 0 $B$J$i(B x_i $B$O(B vGlobal $B$X(B.
%  ii [vlist vGlobal sublist] toGrRing /ii set
%  [ii jj vlist [(partialEcartGlobalVarX) vGlobal]] ecart.isSameIdeal $B$H;H$&(B.
%>
/getGrRing {
  /arg1 set
  [/vlist /vw_vector /ans /vGlobal /sublist /newvlist
   /dlist /tt /i /u /v /k
   ] pushVariables
  [
    /vlist arg1 0 get def
    /vw_vector arg1 1 get def

    vlist isString { [vlist to_records pop] /vlist set } { } ifelse
    vlist { toString } map /vlist set
% dlist $B$O(B [(Dx) (Dy) (Dz)] $B$N%j%9%H(B.
    vlist { /tt set [@@@.Dsymbol tt] cat } map /dlist set  

    /newvlist [ ] def /sublist [ ] def /vGlobal [ ] def
% $B2D49$J?7$7$$JQ?t$r(B newvlist $B$X(B. $BCV49I=$r(B sublist $B$X(B.
    0 1 vlist length 1 sub {
      /i set
%  (u,v) $B$O(B (x_i, Dx_i) $B$KBP$9$k(B weight vector 
      /u vlist i get , vw_vector getGrRing.find  def
      u -1 gt {
        vw_vector , u 1 add , get /u set
      }  { /u 0 def } ifelse

      /v dlist i get , vw_vector getGrRing.find  def
      v -1 gt {
        vw_vector , v 1 add , get /v set
      }  { /v 0 def } ifelse
      u to_int32 /u set , v to_int32 /v set

      u v add , 0  gt {
        newvlist [vlist i get]  join /newvlist set
      } {  } ifelse
      u 0 lt {
        vGlobal [vlist i get] join /vGlobal set
      } {  } ifelse
    } for

    newvlist { /tt set [ [@@@.Dsymbol tt] cat [tt (')] cat ] } map
    /sublist set

    /ans [ vlist , newvlist { /tt set [tt (')] cat } map , join  from_records
           vGlobal sublist] def
    /arg1 ans def
  ] pop
  popVariables
  arg1
} def

%<
% Usages: a uset getGrRing.find index
%>
/getGrRing.find {
   /arg2 set /arg1 set
   [/a /uset /ans /i]  pushVariables
   [
     /a arg1 def /uset arg2 def
     /ans -1 def
     { /ans -1 def
       0 1 , uset length 1 sub {
         /i set
         a tag , uset i get tag eq {
           a , uset i get eq {
             /ans i def  exit
           } { } ifelse 
         } { } ifelse
       } for
       exit
     } loop
     /arg1 ans def
   ] pop
   popVariables
   arg1
} def

%<
% Usages: g1 g2 isSameGrRing bool
%  g1, g2 $B$O(B getGrRing $B$NLa$jCM(B.
%>
/isSameGrRing {
  /arg2 set /arg1 set
  [/g1 /g2 /ans] pushVariables
  [
    /g1 arg1 def /g2 arg2 def
    {
       /ans 1 def
       g1 0 get , g2 0 get eq { } { /ans 0 def exit } ifelse
       exit
       g1 1 get , g2 1 get eq { } { /ans 0 def exit } ifelse
    } loop
    /arg1 ans def
  ] pop 
  popVariables
  arg1
} def

%<
% Usages:  [[ii i_vw_vector] [jj j_vw_vector] vlist] isSameInGrRing_h
% It computes gb.
%>
/isSameInGrRing_h {
  /arg1 set
  [/ii /i_vw_vector /jj /j_vw_vector /vlist 
   /i_gr /j_gr /rrule /ans] pushVariables
  [
    /ii arg1 [0 0] get def
    /i_vw_vector arg1 [0 1] get def
    /jj arg1 [1 0] get def
    /j_vw_vector arg1 [1 1] get def
    /vlist arg1 2 get def
    {
      [vlist i_vw_vector] getGrRing /i_gr set
      [vlist j_vw_vector] getGrRing /j_gr set
      i_gr j_gr isSameGrRing {  } { /ans [0 [i_gr j_gr]] def exit} ifelse

% bug: in case of module
      [i_gr 0 get , ring_of_differential_operators 0] define_ring

% H $B$r(B 1 $B$K(B.
      /rrule [ [@@@.Hsymbol . (1).] ] def

      i_gr 2 get length 0 eq { 
      } {   
        rrule i_gr 2 get  { { . } map } map join /rrule set
      } ifelse
      ii { toString . rrule replace toString } map /ii set
      jj { toString . rrule replace toString } map /jj set

      [ii jj i_gr 0 get , i_gr 1 get] ecartd.isSameIdeal_h /ans set
      [ans [i_gr] rrule ecartd.isSameIdeal_h.failed]  /ans set

      exit
    } loop
    /arg1 ans def
  ] pop 
  popVariables
  arg1 
} def

/test1.isSameInGrRing_h {
  [(parse) (data/test8-data.sm1) pushfile] extension

  cone.gblist 0 get (initial) getNode 2 get /ii set
  cone.gblist 0 get (weight) getNode [2 0 2] get    /iiw set

  cone.gblist 1 get (initial) getNode 2 get /jj set
  cone.gblist 1 get (weight) getNode [2 0 2] get    /jjw set

  (Doing   [ [ii iiw] [jj jjw] cone.vv ] isSameInGrRing_h /ff set) message  
  [ [ii iiw] [jj jjw] cone.vv ] isSameInGrRing_h /ff set

  ff pmat

} def


%<
% Usages: i j isSameCone_h.0  [bool, ...]
% $B%F%9%HJ}K!(B.  (data/test8.sm1) run  (data/test8-data.sm1) run 0 1 isSameCone_h.0
% gb $B$r:FEY7W;;$9$k(B stand alone $BHG(B.  gr(Local ring) $B$GHf3S(B.
%>
/isSameCone_h.0 {
  /arg2 set /arg1 set
  [/i /j /ans /ii /iiw /jj /jjw] pushVariables
  [ 
    /i arg1 def /j arg2 def
    i to_int32 /i set , j to_int32 /j set
    cone.debug { (Comparing ) messagen [i j]  message } { } ifelse

    cone.gblist i get (initial) getNode 2 get /ii set
    cone.gblist i get (weight) getNode [2 0 2] get    /iiw set
 
    cone.gblist j get (initial) getNode 2 get /jj set
    cone.gblist j get (weight) getNode [2 0 2] get    /jjw set

    [ [ii iiw] [jj jjw] cone.vv ] isSameInGrRing_h /ans set

    ans /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
% Usages: [ii vv i_vw_vector] getGbInGrRing_h [ii_gr  i_gr]
% Get Grobner Basis of ii in the graded ring.
% The graded ring is obtained automatically from vv and i_vw_vector.
% ii_gr is the Grobner basis. i_gr is the output of getGrRing.
% cf. isSameInGrRing_h,   ecart.isSameIdeal_h with [(noRecomputation) 1]
%>
/getGbInGrRing_h {
  /arg1 set
  [/ii /i_vw_vector /vlist  /rng /vv /vvGlobal /wv /iigg
   /i_gr  /rrule /ans] pushVariables
  [
    /ii arg1 0 get def
    /vlist arg1 1 get def
    /i_vw_vector arg1 2 get def
    [vlist i_vw_vector] getGrRing /i_gr set

% bug: in case of module
    [i_gr 0 get , ring_of_differential_operators 0] define_ring

% H $B$r(B 1 $B$K(B.
    /rrule [ [@@@.Hsymbol . (1).] ] def

    i_gr 2 get length 0 eq { 
    } {   
      rrule i_gr 2 get  { { . } map } map join /rrule set
    } ifelse
    /vvGlobal i_gr 1 get def
    /vv i_gr 0 get def

    ii { toString . rrule replace toString } map /ii set

    [vv vvGlobal] ecart.stdBlockOrder /wv set
      vvGlobal length 0 eq {
      /rng [vv wv ] def
    }{
      /rng [vv wv [(partialEcartGlobalVarX) vvGlobal]] def
    } ifelse
    /save-cone.autoHomogenize ecart.autoHomogenize def
    /ecart.autoHomogenize 0 def
    [ii] rng join  ecartd.gb  /iigg set
    save-cone.autoHomogenize /ecart.autoHomogenize set
    /ans [iigg 0 get i_gr] def
    /arg1 ans def
  ] pop 
  popVariables
  arg1 
} def

/test1.getGbInGrRing_h {
  [(parse) (data/test8-data.sm1) pushfile] extension

  cone.gblist 0 get (initial) getNode 2 get /ii set
  cone.gblist 0 get (weight) getNode [2 0 2] get    /iiw set
  [ii cone.vv iiw] getGbInGrRing_h /ff1 set

  cone.gblist 1 get (initial) getNode 2 get /jj set
  cone.gblist 1 get (weight) getNode [2 0 2] get    /jjw set
  [jj cone.vv jjw] getGbInGrRing_h /ff2 set

  (ff1 and ff2) message

} def


%<
% setGrGblist
%  cone.grGblist $B$r@_Dj$9$k(B.
%>
/setGrGblist {
  [/ii /ww /gg] pushVariables
  [
    cone.gblist { 
      /gg set
      gg (initial) getNode 2 get /ii set
      gg (weight) getNode [2 0 2] get /ww set
      [ii cone.vv ww] getGbInGrRing_h 
    } map /cone.grGblist set
  ] pop
  popVariables 
} def

%<
% Usages: i j isSameCone_h.2  [bool, ...]
% gb $B$r:FEY7W;;$7$J$$(B.
%>
/isSameCone_h.2 {
  /arg2 set /arg1 set
  [/i /j /ans /ii /iiw /jj /jjw] pushVariables
  [ 
    /i arg1 def /j arg2 def
     i to_int32 /i set , j to_int32 /j set
    (cone.grGblist) boundp { } { setGrGblist } ifelse
    cone.debug { (Comparing ) messagen [i j]  message } { } ifelse

    cone.grGblist i get /ii set
    cone.grGblist j get /jj set

    ii 1 get ,  jj 1 get isSameGrRing {  } 
    { /ans [0 [ii 1 get jj 1 get]] def exit} ifelse

    [ii 0 get , jj 0 get cone.vv [[(noRecomputation) 1]] ] 
    ecartd.isSameIdeal_h /ans set
    [ans [ii 1 get] ii 1 get , ecartd.isSameIdeal_h.failed]  /ans set

    ans /arg1 set
  ] pop
  popVariables
  arg1
} def

%<
%  test1.isSameCone_h.2 $B$O(B cone.grGblist $B$K(B initial $B$N(B gb $B$r(B graded ring
%  $B$G$^$:7W;;$7(B, $B$=$l$+$i(B ideal $B$NHf3S$r$*$3$J$&(B. isSameCone_h.1 $B$KHf$Y$F(B
%  gb $B$N:FEY$N7W;;$,$J$$$N$G7P:QE*(B. 
%> 
/test1.isSameCone_h.2 {
  /cone.loaded boundp { }
  {
    [(parse) (cohom.sm1) pushfile] extension
    [(parse) (dhecart.sm1) pushfile] extension
    /cone.loaded 1 def
  } ifelse
  %[(parse) (cone.sm1) pushfile] extension
  [(parse) (data/test8-data.sm1) pushfile] extension
  setGrGblist  
  (cone.grGblist is set.) message
  0 1 isSameCone_h.2 pmat
} def

%<
% 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
% $B$dBg0hJQ?t$K;H$&(B.
% cone.gblist, cone.fan $B$,@5$7$/@_Dj$5$l$F$$$k$3$H(B.
% (setGrGblist $B$r<B9T:Q$G$"$k$3$H(B. $B<+F0<B9T$5$l$k$,(B... )
%
%>

/isSameCone_h {  isSameCone_h.2 } def

%<
% Usages: genDhcone.init
%   dhcone.checked (dehomogenized $B:Q$N(B cone$BHV9f(B),  dhcone.unchecked $B$N=i4|2=(B.  
%>
/genDhcone.init {
  /dhcone.checked [ ] def
  /dhcone.unchecked [
     0 1 cone.fan length 1 sub {
        to_univNum
     } for
  ] def
} def

%<
% Usages: k genDhcone dhcone
% 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).
%
% $B%F%9%H(B1.  (data/test14.sm1) run (data/test14-data.sm1) run
%          genDhcone.init
%          0 genDhcone /ff set
%>

/genDhcone {
  /arg1 set
  [/k /facets /merged /nextcid /nextfid /coneid
      /newfacets /newmerged /newnextcid /newnextfid /newconeid /vv
   /i /j /p /q /rr /cones /differentC
  ] pushVariables
  [
    /k arg1 def
    /facets [ ] def /merged [ ] def /nextcid [ ] def 
    /nextfid [ ] def /coneid [ ] def
    /cones [ ] def
    /differentC [ ] def

    k to_univNum /k set
    
    {
% Step1. cone.fan[k] $B$r(B $B2C$($k(B.  new... $B$X=i4|%G!<%?$r=q$-9~$`(B.
     cone.debug {(Step 1. Adding ) messagen k messagen (-th cone.) message} { } ifelse
      cones [k to_univNum] join /cones set
      cone.fan k get , (facets) getNode 2 get /vv set
      /newfacets [ ] vv join def

      cone.fan k get , (nextcid) getNode 2 get /vv set
      /newnextcid [ ] vv join def
  
      cone.fan k get , (nextfid) getNode 2 get /vv set
      /newnextfid [ ] vv join def

% newmerged $B$O$^$:(B 0 $B$G$&$a$k(B.  0 : $B$^$@D4$Y$F$J$$(B.
% 1 : merged $B$G>C$($?(B. 2 : boundary. 3 : $B$H$J$j$O0[$J$k(B. 
% [ ] join $B$r$d$C$F(B $B%Y%/%H%k$N(B clone $B$r:n$k(B.
      cone.fan k get , (flipped) getNode 2 get /vv set
      /newmerged [ ] vv join def
      0 1 , newmerged length 1 sub {
         /i set
         newmerged i get , (2).. eq { }
         { newmerged i (0).. put } ifelse
      } for
% newconeid $B$O(B k $B$G$&$a$k(B.
      /newconeid newfacets length newVector { pop k to_univNum } map def

% merged $B$H(B newmerged $B$r(B cone $B$NNY@\4X78$N$_$G99?7$9$k(B.
% $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.
% merged $B$N(B i $BHVL\(B $B$H(B newmerged $B$N(B j $BHVL\$GHf3S(B.
      0 1 , merged length 1 sub {
        /i set
        0 1 , newmerged length 1 sub {
          /j set
          merged i get , (0).. eq , 
          newmerged j get , (0).. eq , and
          nextcid i get , k to_univNum eq , and
          {
             facets i get , newfacets j get , add isZero {
% merged[i], newmerged[j] $B$K(B 1 $B$rF~$l$F>C$9(B.
% $B>e$NH=Dj$O(B nextfid, newnextfid $B$rMQ$$$F$b$h$$$N$G$O(B? 
               merged i (1).. put
               newmerged j (1).. put
             } {  } ifelse
          } { } ifelse
        } for
      } for

% Step2. $B7k9g$7$F$+$i(B, $B$^$@D4$Y$F$J$$(B facet $B$rC5$9(B.
      cone.debug { (Step 2. Joining *** and new***) message } { } ifelse
      /facets facets newfacets join def 
      /merged merged newmerged join def 
      /nextcid nextcid newnextcid join def 
      /nextfid nextfid newnextfid join 
      /coneid  coneid newconeid join def

      cone.debug{ (   Checking facets.) message } { } ifelse
      /k null def
      0 1 , merged length 1 sub {
        /i set
        % i message
        merged i get (0).. eq {
% i $BHVL\$r$^$@D4$Y$F$$$J$$(B.
          coneid i get ,  /p set  
          nextcid i get , /q set
          cone.debug { [p q] message } {  } ifelse
          q (0).. ge {
% cone.fan [p] $B$H(B cone.fan [q] $B$N(B initial $B$rHf3S$9$k(B.
% $BF1$8$J$i(B k $B$r@_Dj(B. exit for. $B0c$($P(B merged[i] = 3 ($B0c$&(B) $B$rBeF~(B.
% differentC $B$O$9$G$K(B $B8=:_$N(B dhcone $B$H0c$&$H(B check $B$5$l$?(B cone $BHV9f(B.
% dhcone.checked $B$O(B dhcone $B$,$9$G$K@8@.$5$l$F$$$k(B cone $BHV9f$N%j%9%H(B.
% $B$3$l$K$O$$$C$F$$$F$b0c$&(B.
            q differentC memberQ , q dhcone.checked memberQ , or 
            { /rr [0 ] def }
            { p q isSameCone_h /rr set } ifelse

            rr 0 get 1 eq {
              cone.debug { (Found next cone. ) message } { } ifelse
              /k q to_univNum def exit
            } {
              cone.debug { ( It is a different cone. ) message } { } ifelse
              differentC [ q ]  join /differentC set
              merged i (3).. put
            } ifelse
          } {  } ifelse
        } {  } ifelse
      } for

      k tag 0 eq { exit } {  } ifelse
   } loop

   [(-1)..] cones join shell rest /cones set
%     dhcone.checked, dhcone.unchecked $B$r99?7(B. 
   dhcone.checked cones join /dhcone.checked set
   dhcone.unchecked cones setMinus /dhcone.unchecked set

   [(dhcone) [ ]
     [
       [(cones) [ ] cones] arrayToTree
       [(facets) [ ] facets] arrayToTree
       [(merged) [ ] merged] arrayToTree
       [(nextcid) [ ] merged] arrayToTree
       [(nextfid) [ ] merged] arrayToTree
       [(coneid) [ ] merged] arrayToTree
     ] 
   ] arrayToTree /arg1 set
  ] pop
  popVariables
  arg1
} def


%<
% Usages: dhCones_h
% cone.fan $B$O(B doubly homogenized (local) $B$G@8@.$5$l$?(B Grobner fan. 
% 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.
%
% $B%F%9%H(B1.  (data/test14.sm1) run (data/test14-data.sm1) run
%          dhCones_h
%          test22
%>
/dhCones_h {
  (cone.grGblist) boundp { } {setGrGblist} ifelse
  genDhcone.init
  /dhcone.fan [ ] def
  {
     (-----------------------------------------) message
     (#dhcone.unchecked = ) messagen dhcone.unchecked length message
     dhcone.unchecked length 0 eq { exit } { } ifelse
     dhcone.fan 
     [ dhcone.unchecked 0 get , genDhcone ] join /dhcone.fan set
     (#dhcone.fan = ) messagen dhcone.fan length message
  } loop
  dhcone.fan
} def

% Todo: print, save functions.  Representative of weight & init.