(*    Mixed bag of Input parametric rational curves
      and surfaces to be polarized  (11/25/2002)                 *)



(* 

  gbicontnet[{cpoly__},p_,q_,r1_,s1_,r2_,s2_] 


  Computes a rectangular control net of bidegree (p, q) 
  from some polar forms, *)
(*  w.r.t. affine frames ((r1, s1), (r2, s2)  *) 

(*  

The input is a list of polynomials, one for each coordinate,
including the denominator (x, y, z, w);
Each polynomial is a list of monomials; 
Each monomial c U^i V^k is a triple {c, i, k}.

Example: Moebius strip, rectangular, degree (6, 1):

poly = {{{2, 0, 0}, {-10, 2, 0}, {-10, 4, 0}, {2, 6, 0}, {2, 1, 1},
         {-12, 3, 1}, {2, 5, 1}},
        {{8, 1, 0}, {-8, 5, 0}, {8, 2, 1}, {-8, 4, 1}},
        {{1, 0, 1}, {1, 2, 1}, {-1, 4, 1}, {-1, 6, 1}},
        {{1, 0, 0}, {3, 2, 0}, {3, 4, 0}, {1, 6, 0}}};

                                                                      *)
(*

 fcontpoly[{cpoly__},m_,r1_,s1_] 

computes the control polygon of a rational curve of degree m   *)
(* specified by  polynomials in cpoly, w.r.t. frame [r1, s1]   *) 
(*  fast method using a recurrence relation *)


(* 


The input is a list of polynomials, one for each coordinate,
including the denominator (x, y, z, w);
Each polynomial is a list of monomials; 
Each monomial c t^i is a triple {c, i}.


Example: seven-leaved  rose  (degree 8):


poly = {{{7, 1}, {-35, 3}, {21, 5}, {-1, 7}},
         {{7, 2}, {-35, 4}, {21, 6}, {-1, 8}},
         {{1, 0}, {4, 2}, {6, 4}, {4, 6}, {1, 8}}}

                                                          *)

(*

 fgcontnet[{cpoly__},m_,reftrig_] 
         
     Computes a triangular control net of degree m  *)
(*  w.r.t. reference triangle                               *)
(*  (r, s, t) = ((r1, r2, r3), (s1, s2, s3), (t1, t2, t3))  *) 
(*  fast method using a recurrence relation *)

(*

The input is a list of polynomials, one for each coordinate,
including the denominator (x, y, z, w);
Each polynomial is a list of monomials; 
Each monomial c U^i V^k is a triple {c, i, k}.

Example, sphere, total degree 2:

sphpoly = {{{2, 1, 0}, {-2, 2, 0}, {-2, 1, 1}},
        {{2, 0, 1}, {-2, 0, 2}, {-2, 1, 1}},
        {{2, 1, 0}, {2, 0, 1}, {-2, 1, 1}, {-1, 0, 0}},
        {{2, 2, 0}, {2, 0, 2}, {-2, 1, 0}, {-2, 0, 1}, {2, 1, 1}, {1, 0, 0}}};


                                                                         *)

(*************************************************************************)


A  = {1, 2, 3, 4, 5, 6, 7, 8};
B  = {3, 5, 7};


S = {1, 2, 3, 4};

S = {1, 2, 3, 4, 5};

S = {1, 2, 3, 4, 5, 6, 7, 8};

zozo = subsets[S, 1]
zozo = subsets[S, 2]


SS = {{2, 3}, {2, 4}, {3, 4}};
a = 1;

sss = {1, 2, 3, 4};
sss = {1, 2, 3, 4, 5, 6};

mon1 = {1, 2, 2};
mon2 = {1, 4, 0};
mon3 = {1, 3, 1};
mon4 = {1, 0, 4};
mon5 = {1, 1, 1};
mon6 = {1, 1, 2};
mon7 = {1, 0, 0};
mon8 = {1, 2, 4};
mon9 = {1, 1, 4};
mon10 = {1, 1, 3};
mon11 = {1, 1, 1};
mon12 = {1, 1, 4};
mon13 = {1, 3, 0};


zog = bipolarmon[mon1, 2, 2]
zog = bipolarlismon[mon1, 3, 3]
zog = bipolarmon[mon13, 4, 0]

zog = polarmon[mon2,  4]

zog = polarlismon[mon2,  4]

poly1 = {{1, 1, 0}, {-1/3, 3, 0}, {1, 1, 2}};
poly2 = {{1, 0, 1}, {-1/3, 0, 3}, {1, 2, 1}};
poly3 = {{1, 2, 0},  {-1, 0, 2}};

zig = polarize[poly1,3]
zig = polarize[poly2,3]
zig = polarize[poly3,3]
zig = polarizelis[poly1,3]
zig = bipolarize[poly1,3,2]
zig = bipolarize[poly2,2,3]
zig = bipolarize[poly3,2,2]
zig = bipolarizelis[poly1,3,2]

val1 = evalpoly[aff1, U, V]

(*  Enneper surface   *)

poly = {{{1, 1, 0}, {-1/3, 3, 0}, {1, 1, 2}}, 
        {{1, 0, 1}, {-1/3, 0, 3}, {1, 2, 1}},
        {{1, 2, 0},  {-1, 0, 2}}};

polarform = polarizesurf[poly,3]

bipolarform = bipolarizesurf[poly,3, 3]


recnet =  bicontnet[poly, 3, 3]

recnet =   {{0, 0, 0, 1}, {0, 1/3, 0, 1}, {0, 2/3, -1/3, 1}, {0, 2/3, -1, 1},
   {1/3, 0, 0, 1}, {1/3, 1/3, 0, 1}, {2/3, 2/3, -1/3, 1}, {2/3, 0, 1/3, 1},
   {2/3, 2/3, 1/3, 1}, {2/3, 0, 1, 1}}


(*  poly for Mobius stip  )

poly = {{{2, 0, 0}, {-10, 2, 0}, {-10, 4, 0}, {2, 6, 0}, {2, 1, 1},
         {-12, 3, 1}, {2, 5, 1}},
        {{8, 1, 0}, {-8, 5, 0}, {8, 2, 1}, {-8, 4, 1}},
        {{1, 0, 1}, {1, 2, 1}, {-1, 4, 1}, {-1, 6, 1}},
        {{1, 0, 0}, {3, 2, 0}, {3, 4, 0}, {1, 6, 0}}};

poly1 = {{2, 0, 0}, {-10, 2, 0}, {-10, 4, 0}, {2, 6, 0}, {2, 1, 1},
         {-12, 3, 1}, {2, 5, 1}}

poly2 =  {{8, 1, 0}, {-8, 5, 0}, {8, 2, 1}, {-8, 4, 1}}

poly4 =  {{1, 0, 0}, {3, 2, 0}, {3, 4, 0}, {1, 6, 0}}
poly5 =  {{1, 0, 0}, {1, 6, 0}}


aff1 = bipolarizelis[poly1, 6, 1]

aff2 = bipolarizelis[poly2, 6, 1]
aff4 = bipolarizelis[poly4, 6, 1]
aff5 = bipolarizelis[poly5, 6, 1]



U = {-1, 1, 1, 1, 1, 1};
V = {1};

valo = evalpoly[aff1, U, V]
valo = evalpoly[aff2, U, V]
valo = evalpoly[aff4, U, V]
valo = evalpoly[aff5, U, V]


poly0 = {{1, 0, 0}}
net = contnet[poly, 3]

(* Steiner Roman surface  *)

poly = {{{2, 0, 1}}, {{2, 1, 0}}, {{2, 1, 1}}, 
        {{1, 2, 0}, {1, 0, 2}, {1, 0, 0}}};

net = contnet[poly, 2]

pol1 = {{2, 1, 1}};
zag = polarizelis[pol1, 2]
zap = evalpoly[zag, U, V]

U = {1, 0};
V = {0, 1};


(*  Monkey saddle   *)

poly = {{{1, 1, 0}}, {{1, 0, 1}},
        {{1, 3, 0}, {-3, 1, 2}}};
 
net = contnet[poly, 3]

net = bicontnet[poly, 3, 2]

monknet =  {{0, 0, 0, 1}, {0, 1/2, 0, 1}, {0, 1, 0, 1}, {1/3, 0, 0, 1},
   {1/3, 1/2, 0, 1}, {1/3, 1, -1, 1}, {2/3, 0, 0, 1}, {2/3, 1/2, 0, 1},
   {2/3, 1, -2, 1}, {1, 0, 1, 1}, {1, 1/2, 1, 1}, {1, 1, -2, 1}}



(*  Whitney's umbrella   *)

poly = {{{2, 1, 1}}, {{1, 1, 0}}, {{1, 0, 2}}};

net = contnet[poly, 2]


(* cross-cap surface  (image of projective plane)        *)

poly = {{{2, 2, 1}, {2, 0, 3}, {-2, 0, 1}}, 
        {{8, 1, 1}}, 
        {{4, 2, 0}, {-4, 0, 2}},  
        {{1, 4, 0}, {1, 0, 4}, {2, 2, 2}, {2, 2, 0}, {2, 0, 2}, {1, 0, 0}}};

net = contnet[poly, 4]

net =   {{0, 0, 0, 1}, {-1/2, 0, 0, 1}, {-3/4, 0, -1/2, 4/3}, {-1/2, 0, -1, 2}, 
   {0, 0, -1, 4}, {0, 0, 0, 1}, {-1/2, 2/3, 0, 1}, {-3/4, 1, -1/2, 4/3}, 
   {-1/2, 1, -1, 2}, {0, 0, 1/2, 4/3}, {-1/4, 1, 1/2, 4/3}, 
   {-1/3, 4/3, 0, 2}, {0, 0, 1, 2}, {0, 1, 1, 2}, {0, 0, 1, 4}}


(* algebraic cross-cap surface, second parameterization   *)  

poly = {{{2, 0, 1}, {-2, 4, 1}}, 
        {{4, 1, 1}, {4, 3, 1}}, 
        {{2, 4, 0}, {4, 2, 0}, {2, 0, 0}},  
        {{1, 4, 2}, {1, 4, 0}, {1, 0, 2}, {1, 0, 0}}};

net = contnet[poly, 6]

(*  Plucker's conoid     *)

poly = {{{1, 0, 1}, {-1, 4, 1}},
        {{2, 1, 1}, {2, 3, 1}},
        {{4, 1, 0}, {-4, 3, 0}},
        {{1, 4, 0}, {2, 2, 0}, {1, 0, 0}}};

net = contnet[poly, 5]


(*  Right  conoid  degree 3   *)


poly = {{{1, 0, 1}, {-1, 2, 1}},
        {{2, 1, 1}},
        {{2, 1, 0}},
        {{1, 2, 0}, {1, 0, 0}}};


net = contnet[poly, 3]




(* Torus , degree 4                 *)

poly  =  {{{-aa, 2, 2}, {2 * bb, 2, 1}, {-aa, 2, 0}, {aa, 0, 2}, 
           {-2 * bb, 0, 1}, {aa, 0, 0}},
          {{2 * aa, 1, 2}, {-4 * bb, 1, 1}, {2 * aa, 1, 0}},
          {{-cc, 2, 2}, {-cc, 0, 2}, {cc, 2, 0}, {cc, 0, 0}},
          {{1, 2, 2}, {1, 2, 0}, {1, 0, 2}, {1, 0, 0}}};

 reftrig1 = {{-1, 1, 1}, {-1, -1, 3}, {1, 1, -1}};

net = gcontnet[poly, 4, reftrig1]; tornets2 = InputForm[net];

tornets2 =  {{0, (4*aa - 4*bb)/4, 0, 4}, {2*aa - 2*bb, 0, 2*cc, 0}, {0, bb, 0, 4/3},
   {2*aa + 2*bb, 0, 2*cc, 0}, {0, (-4*aa - 4*bb)/4, 0, 4},
   {(2*aa - 2*bb)/2, (2*aa - 2*bb)/2, 0, 2},
   {(3*((4*aa)/3 - (4*bb)/3))/2, (3*((-2*aa)/3 + (2*bb)/3))/2, cc, 2/3},
   {(3*((2*aa)/3 + (2*bb)/3))/2, (3*((-2*aa)/3 + (2*bb)/3))/2, 2*cc, 2/3},
   {0, (-2*aa - 2*bb)/2, cc, 2}, {(3*((8*aa)/3 - (8*bb)/3))/4, 0, 0, 4/3},
   {(3*((2*aa)/3 - (2*bb)/3))/4, (3*((-4*aa)/3 + (4*bb)/3))/4, cc/2, 4/3},
   {0, -aa, 2*cc, 4/3}, {(2*aa - 2*bb)/2, (-2*aa + 2*bb)/2, 0, 2},
   {0, (-2*aa + 2*bb)/2, cc, 2}, {0, (-4*aa + 4*bb)/4, 0, 4}}

(*   aa = 2,  bb = 1  cc = 1    *)


tornets2 =  {{0, 1, 0, 4}, {2, 0, 2, 0}, {0, 1, 0, 4/3}, {6, 0, 2, 0}, {0, -3, 0, 4},
   {1, 1, 0, 2}, {2, -1, 1, 2/3}, {3, -1, 2, 2/3}, {0, -3, 1, 2},
   {2, 0, 0, 4/3}, {1/2, -1, 1/2, 4/3}, {0, -2, 2, 4/3}, {1, -1, 0, 2},
   {0, -1, 1, 2}, {0, -1, 0, 4}}

tornets3 = maptohat[tornets2]

tornets3 =  {{0, 4, 0, 4}, {2, 0, 2, 0}, {0, 4/3, 0, 4/3}, {6, 0, 2, 0},
   {0, -12, 0, 4}, {2, 2, 0, 2}, {4/3, -2/3, 2/3, 2/3}, {2, -2/3, 4/3, 2/3},
   {0, -6, 2, 2}, {8/3, 0, 0, 4/3}, {2/3, -4/3, 2/3, 4/3},
   {0, -8/3, 8/3, 4/3}, {2, -2, 0, 2}, {0, -2, 2, 2}, {0, -4, 0, 4}}


(* polynomials for torus after change of var u -> v/u,  v -> 1/u   *)


 pol1 = (aa*(u^2 + 1) - 2*bb*u)(u^2 - v^2);
 pol2 = 2*u*v*(aa*(u^2 + 1) - 2*bb*u);
 pol3 = cc*(u^2 -1)*(u^2 + v^2);
 pol4 = (u^2 + 1)*(u^2 + v^2);

 p1 = Expand[pol1]; pp1 = InputForm[p1];
 p2 = Expand[pol2]; pp2 = InputForm[p2];
 p3 = Expand[pol3]; pp3 = InputForm[p3];
 p4 = Expand[pol4]; pp4 = InputForm[p4];

 pp1 =   aa*u^2 - 2*bb*u^3 + aa*u^4 - aa*v^2 + 2*bb*u*v^2 - aa*u^2*v^2
 pp2 = 2*aa*u*v - 4*bb*u^2*v + 2*aa*u^3*v
 pp3 =  -(cc*u^2) + cc*u^4 - cc*v^2 + cc*u^2*v^2
 pp4 = u^2 + u^4 + v^2 + u^2*v^2

 poly1a = {{{aa, 2, 0},  {-2*bb, 3, 0},  {aa, 4, 0},
            {-aa, 0, 2},  {2*bb, 1, 2},  {-aa, 2, 2}},
           {{2*aa, 1, 1},  {-4*bb, 2, 1}, {2*aa, 3, 1}},
           {{-cc, 2, 0},  {cc, 4, 0}, {-cc, 0, 2},  {cc, 2, 2}},
           {{1, 2, 0}, {1, 4, 0}, {1, 0, 2}, {1, 2, 2}}}

 ynet = gcontnet[poly1a, 4, reftrig1]; tornetsa = InputForm[ynet];

 theta1 =   {{0, 1, 0, 4}, {0, 1, 1, -2}, {0, 2, 2, 4/3}, {0, 3, 1, -2}, {0, 3, 0, 4},
   {-2, 0, -2, 0}, {1, 3, 2, 2/3}, {3/2, 3, 1/2, -4/3}, {3, 3, 0, 2},
   {0, 1, 0, 4/3}, {6, 3, 1, -2/3}, {6, 0, 0, 4/3}, {-6, 0, -2, 0},
   {3, -3, 0, 2}, {0, -3, 0, 4}}

 theta1 = maptohat[theta1]

 theta1 =   {{0, 4, 0, 4}, {0, -2, -2, -2}, {0, 8/3, 8/3, 4/3}, {0, -6, -2, -2},
   {0, 12, 0, 4}, {-2, 0, -2, 0}, {2/3, 2, 4/3, 2/3}, {-2, -4, -2/3, -4/3},
   {6, 6, 0, 2}, {0, 4/3, 0, 4/3}, {-4, -2, -2/3, -2/3}, {8, 0, 0, 4/3},
   {-6, 0, -2, 0}, {6, -6, 0, 2}, {0, -12, 0, 4}}

  triga = {{1, 0, 0}, {0, 0, 1}, {0, 1/2, 1/2}};
  
  tnet1 = newcnet3[theta1, 4, triga]
  tnet3 = InputForm[tnet1]
  tnet3 = 
  {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 2/3, -1/3, 1/3}, {0, 1, -1, 1},
   {0, 4, 0, 4}, {0, 0, 0, 0}, {-2/3, 0, 0, 0}, {-2/3, 1, -1/3, 1/3},
   {-2, 0, -2, 0}, {0, -2/3, -1/3, 1/3}, {-2, -1/3, -1/3, 1/3},
   {0, 4/3, 0, 4/3}, {0, -3, -1, 1}, {-6, 0, -2, 0}, {0, -12, 0, 4}}


 trigy = {{0, 0, 1}, {1, 1, -1}, {-1, 1, 1}};

 ynet = gcontnet[poly1a, 4, trigy]; tornetsa = InputForm[ynet];

 tornetsa =  {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, -aa, -cc, 1/3}, {0, -aa - bb, -cc, 1},
   {0, (-4*aa - 4*bb)/4, 0, 4}, {0, 0, 0, 0}, {-aa/3, 0, 0, 0},
   {3*((-2*aa)/3 - (2*bb)/3), 3*(-aa/3 + bb/3), -cc, 1/3},
   {-2*aa - 2*bb, 0, -2*cc, 0}, {0, aa, -cc, 1/3},
   {3*((-2*aa)/3 + (2*bb)/3), 3*(aa/3 + bb/3), -cc, 1/3}, {0, bb, 0, 4/3},
   {0, aa - bb, -cc, 1}, {-2*aa + 2*bb, 0, -2*cc, 0},
   {0, (4*aa - 4*bb)/4, 0, 4}}


 tornetsa = {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, -2, -1, 1/3}, {0, -3, -1, 1},
   {0, -3, 0, 4}, {0, 0, 0, 0}, {-2/3, 0, 0, 0}, {-6, -1, -1, 1/3},
   {-6, 0, -2, 0}, {0, 2, -1, 1/3}, {-2, 3, -1, 1/3}, {0, 1, 0, 4/3},
   {0, 1, -1, 1}, {-2, 0, -2, 0}, {0, 1, 0, 4}}



 tornesta2 = maptohat[tornetsa]
 tornetsa3 = InputForm[tornetsa2]

 tornetsa3 = {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, -2/3, -1/3, 1/3}, {0, -3, -1, 1},
   {0, -12, 0, 4}, {0, 0, 0, 0}, {-2/3, 0, 0, 0}, {-2, -1/3, -1/3, 1/3},
   {-6, 0, -2, 0}, {0, 2/3, -1/3, 1/3}, {-2/3, 1, -1/3, 1/3},
   {0, 4/3, 0, 4/3}, {0, 1, -1, 1}, {-2, 0, -2, 0}, {0, 4, 0, 4}}



 trigz = {{0, 0, 1}, {-1, 1, 1}, {-1, -1, 3}};

 znet = gcontnet[poly1a, 4, trigz]; tornetsb = InputForm[znet];

 tornetsb = {{0, 3, 0, 4}, {3, 3, 0, 2}, {6, 0, 0, 4/3}, {3, -3, 0, 2}, {0, -3, 0, 4},
   {0, 3, -1, 1}, {6, 3, -1, 1/3}, {6, -3, -1, 1/3}, {0, -3, -1, 1},
   {0, 2, -1, 1/3}, {2/3, 0, 0, 0}, {0, -2, -1, 1/3}, {0, 0, 0, 0},
   {0, 0, 0, 0}, {0, 0, 0, 0}}
 
 tornetsb2 = maptohat[tornetsb]
 tornetsb3 = InputForm[tornetsb2]

 tornetsb2 =  {{0, 12, 0, 4}, {6, 6, 0, 2}, {8, 0, 0, 4/3}, {6, -6, 0, 2},
   {0, -12, 0, 4}, {0, 3, -1, 1}, {2, 1, -1/3, 1/3}, {2, -1, -1/3, 1/3},
   {0, -3, -1, 1}, {0, 2/3, -1/3, 1/3}, {2/3, 0, 0, 0}, {0, -2/3, -1/3, 1/3},
   {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}}



(* polarizing polynomials in one variable  *)

(* Lemniscate of Bernoulli  *)


mon1 = {1, 4}
mon2 = {1, 1}
mon3 = {1, 0}

poly1 = {{1, 1}, {1, 3}}
poly2 = {{1, 1}, {-1, 3}}
poly3 = {{1, 0}, {1, 4}}


aff1 = cpolarmon[mon1, 4]
aff1 = cpolarlismon[mon1, 4]
aff1 = cpolarmon[mon2, 4]
aff1 = cpolarlismon[mon2, 4]
aff1 = cpolarmon[mon3, 4]
aff1 = cpolarlismon[mon3, 4]

aff1 = cpolarize[poly1, 4]
aff1 = cpolarizelis[poly1, 4]
aff1 = cpolarize[poly2, 4]
aff1 = cpolarizelis[poly2, 4]
aff1 = cpolarize[poly3, 4]
aff1 = cpolarizelis[poly3, 4]


poly = {{{1, 1}, {1, 3}},
        {{1, 1}, {-1, 3}},
        {{1, 0}, {1, 4}}}

cpoly = gcontpoly[poly, 4, 0, 1]

cpoly =   {{0, 0, 1}, {1/4, 1/4, 1}, {1/2, 1/2, 1}, {1, 1/2, 1}, {1, 0, 2}}


(* Folium of Descartes  *)

poly = {{{3, 1}},
        {{3, 2}},
        {{1, 0}, {1, 3}}}

cpoly = gcontpoly[poly, 3, 0, 1]


(* four-leaved rose WRONG! *)

p1 = 4*t*(1 - t^2)^2;
p2 = 8*t^2*(1 - t^2)^2;
p3 = (1 + t^2)^3;
pp1 = Expand[p1]; pol1 = InputForm[pp1]
pp2 = Expand[p2]; pol2 = InputForm[pp2]
pp3 = Expand[p3]; pol3 = InputForm[pp3]

pol1 = 4*t - 8*t^3 + 4*t^5
pol2 =  8*t^2 - 16*t^4 + 8*t^6
pol3 = 1 + 3*t^2 + 3*t^4 + t^6

poly = {{{4, 1}, {-8, 3}, {4, 5}},
        {{8, 2}, {-16, 4}, {8, 6}},
        {{1, 0}, {3, 2}, {3, 4}, {1, 6}}}


cpoly = gcontpoly[poly, 6, 0, 1]

cpoly =   {{0, 0, 1}, {2/3, 0, 1}, {10/9, 4/9, 6/5}, {1, 1, 8/5}, {4/9, 8/9, 12/5},
   {0, 0, 4}, {0, 0, 8}}



(* four-leaved rose  *)

p1 = 4*t*(1 - t^2)^2;
p2 = 8*t^2*(1 - t^2);
p3 = (1 + t^2)^3;
pp1 = Expand[p1]; pol1 = InputForm[pp1]
pp2 = Expand[p2]; pol2 = InputForm[pp2]
pp3 = Expand[p3]; pol3 = InputForm[pp3]

pol1 = 4*t - 8*t^3 + 4*t^5
pol2 = 8*t^2 - 8*t^4
pol3 = 1 + 3*t^2 + 3*t^4 + t^6

poly = {{{4, 1}, {-8, 3}, {4, 5}},
        {{8, 2}, {-8, 4}},
        {{1, 0}, {3, 2}, {3, 4}, {1, 6}}}


cpoly = gcontpoly[poly, 6, 0, 1]

cpoly =    {{0, 0, 1}, {2/3, 0, 1}, {10/9, 4/9, 6/5}, {1, 1, 8/5}, {4/9, 10/9, 12/5},
   {0, 2/3, 4}, {0, 0, 8}}

cpoly = fcontpoly[poly, 6, -1, 1]
 
(* five-leaved  rose   *)

p1 = t*(5 - 10*t^2 + t^4)
p2 = t^2*(5 - 10*t^2 + t^4)
p3 = (1 + t^2)^3

pol1 = 5*t - 10*t^3 + t^5
pol2 = 5*t^2 - 10*t^4 + t^6
pol3 = 1 + 3*t^2 + 3*t^4 + t^6

poly = {{{5, 1}, {-10, 3}, {1, 5}},
        {{5, 2}, {-10, 4}, {1, 6}},
        {{1, 0}, {3, 2}, {3, 4}, {1, 6}}}        

cpoly = gcontpoly[poly, 6, 0, 1]


cpoly =   {{0, 0, 1}, {5/6, 0, 1}, {25/18, 5/18, 6/5}, {5/4, 5/8, 8/5},
   {5/9, 5/9, 12/5}, {-1/6, 0, 4}, {-1/2, -1/2, 8}}


(* seven-leaved  rose   *)

p3 = (1 + t^2)^4
p3 = Expand[p3]; pol3 = InputForm[p3]

pol1 = 7*t - 35*t^3 + 21*t^5 - t^7
pol2 = 7*t^2 - 35*t^4 + 21*t^6 - t^8
pol3 = 1 + 4*t^2 + 6*t^4 + 4*t^6 + t^8


poly = {{{7, 1}, {-35, 3}, {21, 5}, {-1, 7}},
         {{7, 2}, {-35, 4}, {21, 6}, {-1, 8}},
         {{1, 0}, {4, 2}, {6, 4}, {4, 6}, {1, 8}}}

ppol =   {{1, 0}, {4, 2}, {6, 4}, {4, 6}, {1, 8}}

aff1 = cpolarize[ppol, 8]
aff1 = cpolarizelis[ppol, 8]

U = {1, 1, 1, 1, 1, 1, 1, 1}

vall = gevalpoly1[aff1, U]

cpoly = gcontpoly[poly, 8, 0, 1]
cpoly1 = fcontpoly[poly, 8, 0, 1]


cpoly =    {{0, 0, 1}, {7/8, 0, 1}, {49/32, 7/32, 8/7}, {7/5, 21/40, 10/7},
   {35/68, 35/68, 68/35}, {-21/40, 0, 20/7}, {-35/32, -21/32, 32/7},
   {-1, -7/8, 8}, {-1/2, -1/2, 16}}



cpoly2 = gcontpoly[poly, 8, -1, 1]

(* fast method *)
cpoly2 = fcontpoly[poly, 8, -1, 1];

toto = InputForm[cpoly2]

cpoly2 = {{1/2, -1/2, 16}, {8, -6, 0}, {-7/4, 7/2, 16/7}, {-8, 2, 0}, 
   {0, -35/6, 48/35}, {8, 2, 0}, {7/4, 7/2, 16/7}, {-8, -6, 0}, 
   {-1/2, -1/2, 16}}


cpoly3 = fcontpoly[poly, 8, -2, 1]

toto = InputForm[cpoly3]


cpoly3 = {{-278/625, 556/625, 625}, {-233/1000, 13/10, -125}, 
   {2233/1600, -2933/1600, 400/7}, {28/85, -3031/680, -170/7}, 
   {-1225/548, -805/548, 548/35}, {-77/136, 91/34, -68/7}, 
   {133/64, 175/64, 64/7}, {1, 5/8, -8}, {-1/2, -1/2, 16}}




(* three leaved rose    *)

pol1 = 3*t - t^3
pol2 = 3*t^2 - t^4
pol4 = 1 + 2*t^2 + t^4

poly = {{{3, 1}, {-1, 3}},
        {{3, 2}, {-1, 4}},
        {{1, 0}, {2, 2}, {1, 4}}}


cpoly = gcontpoly[poly, 4, 0, 1]

cpoly =  {{0, 0, 1}, {3/4, 0, 1}, {9/8, 3/8, 4/3}, {1, 3/4, 2}, {1/2, 1/2, 4}}


(* Viviani Window    *)


pol1 = 2*t - 2*t^3
pol2 = 4*t^2
pol3 = 1 - t^4
pol4 = 1 + 2*t^2 + t^4

poly = {{{2, 1}, {-2, 3}},
        {{4, 2}},
        {{1, 0}, {-1, 4}},
        {{1, 0}, {2, 2}, {1, 4}}}


cpoly = gcontpoly[poly, 4, 0, 1]


cpoly =   {{0, 0, 1, 1}, {1/2, 0, 1, 1}, {3/4, 1/2, 3/4, 4/3}, {1/2, 1, 1/2, 2},
   {0, 1, 0, 4}}


cpoly = fcontpoly[poly, 4, 0, 1]


(* Lissajous type curve     *)

p1 = (1 - t^2)*(1 - 14*t^2 + t^4);
p2 = 4*t*(1 - t^2)*(1 + t^2);
p3 = (1 + t^2)^3;


pol1 = Expand[p1]; pp1 = InputForm[pol1]
pol2 = Expand[p2]; pp2 = InputForm[pol2]
pol3 = Expand[p3]; pp3 = InputForm[pol3]


pp1 = 1 - 15*t^2 + 15*t^4 - t^6
pp2 =  4*t - 4*t^5
pp3 = 1 + 3*t^2 + 3*t^4 + t^6

poly = {{{1, 0}, {-15, 2}, {15, 4}, {-1, 6}},
        {{4, 1}, {-4, 5}},
        {{1, 0}, {3, 2}, {3, 4}, {1, 6}}}


cpoly = gcontpoly[poly, 6, 0, 1]

cpoly =   {{1, 0, 1}, {1, 2/3, 1}, {0, 10/9, 6/5}, {-5/4, 5/4, 8/5},
   {-5/3, 10/9, 12/5}, {-1, 2/3, 4}, {0, 0, 8}}


(* sphere  *)

sphpoly = {{{2, 1, 0}, {-2, 2, 0}, {-2, 1, 1}},
        {{2, 0, 1}, {-2, 0, 2}, {-2, 1, 1}},
        {{2, 1, 0}, {2, 0, 1}, {-2, 1, 1}, {-1, 0, 0}},
        {{2, 2, 0}, {2, 0, 2}, {-2, 1, 0}, {-2, 0, 1}, {2, 1, 1}, {1, 0, 0}}};


reftrig = {{1, 0, 0}, {0, 1, 0}, {0, 0, 1}};

 cnet = fgcontnet[sphpoly, 2, reftrig];  spnet = InputForm[cnet];

 spnet =   {{0, 0, -1, 1}, {0, 1, 0, 0}, {0, 0, 1, 1}, {1, 0, 0, 0}, {0, 0, 0, 0},
   {0, 0, 1, 1}}



(* Steiner surface  *)
(* 1st quarter  *)

(* x = 2(v - v^2 - uv) = 2v(1 - u - v)
   y = 2(u - u^2 - uv) = 2u(1 - u - v)
   z = 2uv
   w = 2u^2 + 2v^2 -2u -2v +2uv + 1 = u^2 + v^2 + (1 - u - v)^2
*)

steinpoly1 = {{{2, 0, 1}, {-2, 0, 2}, {-2, 1, 1}},
        {{2, 1, 0}, {-2, 2, 0}, {-2, 1, 1}},
        {{2, 1, 1}},
        {{2, 2, 0}, {2, 0, 2}, {-2, 1, 0}, {-2, 0, 1}, {2, 1, 1}, {1, 0, 0}}};


 cnet1 = fgcontnet[steinpoly1, 2, reftrig];  stneth1 = InputForm[cnet1];

 stneth1 =   {{0, 0, 0, 1}, {1, 0, 0, 0}, {0, 0, 0, 1}, {0, 1, 0, 0}, {0, 0, 1, 0},
   {0, 0, 0, 1}};

(* 2nd quarter  *)

(* x = 2(v - v^2 - uv)
   y = 2(-u + u^2 + uv)  
   z = -2uv
   w = 2u^2 + 2v^2 -2u -2v +2uv + 1
*)

steinpoly2 = {{{2, 0, 1}, {-2, 0, 2}, {-2, 1, 1}},
        {{-2, 1, 0}, {2, 2, 0}, {2, 1, 1}},
        {{-2, 1, 1}},
        {{2, 2, 0}, {2, 0, 2}, {-2, 1, 0}, {-2, 0, 1}, {2, 1, 1}, {1, 0, 0}}};


 cnet2 = fgcontnet[steinpoly2, 2, reftrig];  stneth2 = InputForm[cnet2];


 stneth2 =   {{0, 0, 0, 1}, {1, 0, 0, 0}, {0, 0, 0, 1}, {0, -1, 0, 0}, {0, 0, -1, 0},
   {0, 0, 0, 1}}


(* 3rd quarter  *)

(* x = 2(-v + v^2 + uv)
   y = 2(u - u^2 - uv)  
   z = -2uv
   w = 2u^2 + 2v^2 +2u -2v -2uv + 1
*)

steinpoly3 = {{{-2, 0, 1}, {2, 0, 2}, {2, 1, 1}},
        {{2, 1, 0}, {-2, 2, 0}, {-2, 1, 1}},
        {{-2, 1, 1}},
        {{2, 2, 0}, {2, 0, 2}, {-2, 1, 0}, {-2, 0, 1}, {2, 1, 1}, {1, 0, 0}}};

 cnet3 = fgcontnet[steinpoly3, 2, reftrig];  stneth3 = InputForm[cnet3];

 stneth3 =  {{0, 0, 0, 1}, {-1, 0, 0, 0}, {0, 0, 0, 1}, {0, 1, 0, 0}, {0, 0, -1, 0},
   {0, 0, 0, 1}}


(* 4th quarter  *)

(* x = 2(-v + v^2 + uv)
   y = 2(-u + u^2 + uv)  
   z = 2uv
   w = 2u^2 + 2v^2 -2u +2v -2uv + 1
*)


steinpoly4 = {{{-2, 0, 1}, {2, 0, 2}, {2, 1, 1}},
        {{-2, 1, 0}, {2, 2, 0}, {2, 1, 1}},
        {{2, 1, 1}},
        {{2, 2, 0}, {2, 0, 2}, {-2, 1, 0}, {-2, 0, 1}, {2, 1, 1}, {1, 0, 0}}};

 cnet4 = fgcontnet[steinpoly4, 2, reftrig];  stneth4 = InputForm[cnet4];

 stneth4 =    {{0, 0, 0, 1}, {-1, 0, 0, 0}, {0, 0, 0, 1}, {0, -1, 0, 0}, {0, 0, 1, 0},
   {0, 0, 0, 1}}