legacy/oldlocal/ch_local/local_fan.ml

   1 (*
   2 Created Nov 26, 2009
   3 THALES
   4 Compute irreducible cyclic fans.
   5 *)
   6
   7 (*
   8 Needs.
   9 k nodes, k edges {i,i+1},
  10 cyclic map on both
  11 Each node is 2, 2h0 or free.
  12 Each node is flat or not flat.
  13
  14 Each edge is in G or not,
  15 Each edge is 2, 2h0, or free.
  16
  17 *)
  18 open List;;
  19
  20 (* from lpproc.ml *)
  21 let upto =
  22   let rec rangeA a i j  = if (i >= j) then a
  23    else rangeA ((j-1)::a) i (j-1)  in
  24   rangeA [] 0;;
  25
  26
  27 (* minimal fan definition *)
  28
  29 type node_t = N2 | N2h0 | Nfree;;
  30 type edge_t = E2 | E2h0 | Gset;;
  31
  32
  33 type minimal_fan =  {
  34   node : node_t list;
  35   nodeflats : bool list;
  36   edge : edge_t list;
  37 };;
  38
  39 let mk_minimal_fan nhts nflts e =
  40    {
  41    node = nhts;
  42    nodeflats = nflts;
  43    edge = e;
  44    };;
  45
  46
  47 (* generating minimal_fan *)
  48
  49 let node_of_int i = if (i=0) then N2 else if (i=1) then N2h0 else Nfree;;
  50 let edge_of_int i = if (i=0) then E2 else if (i=1) then E2h0 else Gset;;
  51 let bool_of_int i = if (i=0) then false else true;;
  52
  53 let base modulus len k =
  54   let rec baselist modulus len k acc =
  55     if len <= 0 then acc else
  56     if (k=0) then baselist modulus (len-1) 0 (0::acc)
  57     else
  58       baselist modulus (len - 1) (k/modulus) ((k mod modulus) :: acc) in
  59   baselist modulus len k [];;
  60
  61 let rec pow base exp =
  62    if exp = 0 then 1 else base*(pow base (exp-1));;
  63
  64 let mk_x len (a,b,c) (r,s,t) =
  65   let (n,nf,e) = (base a len r, base b len s, base c len t) in
  66   mk_minimal_fan
  67      (map node_of_int n)
  68      (map bool_of_int nf)
  69      (map edge_of_int e);;
  70
  71 (* reading data from a record *)
  72
  73 let kn mf = length (mf.node);;
  74 let sn mf = List.length (filter ((=) Gset) mf.edge);;
  75 let rn mf = kn mf - sn mf;;
  76
  77 let posmod i m =
  78   let j = i mod m in
  79   if (j<0) then j+m else j;;
  80
  81 let part xs i = nth xs (posmod i (length xs));;
  82
  83 let number mf = upto(length mf.edge);;
  84
  85 let g_edge mf i = (part mf.edge i = Gset);;
  86 (* let nong_edge mf i = not(g_edge mf i);; *)
  87 let gminimal_edge mf i = not(part mf.edge i = E2h0);;
  88
  89 let flat_node mf i = part mf.nodeflats i;;
  90 let nonflat_node mf i = not(flat_node mf i);;
  91 let bound_node mf i = not(part mf.node i = Nfree);;
  92
  93 (* irreducibility *)
  94
  95 (* extreme_edge is built into construction of edge types *)
  96
  97 let card mf = (sn mf <= 3) && (3 <= sn mf + rn mf) && (rn mf + 2 * sn mf <= 6);;
  98
  99 let extreme_edge mf = true;;
 100
 101 let flat_exists mf =
 102   let has i = (kn mf <= 4) or
 103        flat_node mf i or flat_node mf (i+1) or flat_node mf (i+2) or flat_node mf (i+3) in
 104   for_all has (number mf);;
 105
 106 let no_triple_flat mf =
 107    let triple_flat i = flat_node mf i && flat_node mf (i+1) && flat_node mf (i+2) in
 108    not (exists triple_flat (number mf));;
 109
 110 let balance mf =
 111    let es = mf.edge in
 112    let has_balance i = (part es i = part es (i+1)) or (part es i = Gset) or (part es (i+1) = Gset) in
 113    for_all has_balance (number mf);;
 114
 115 let g_flat mf =
 116    let gg = g_edge mf in
 117   let nf = nonflat_node mf in
 118   let has_gflat i = nf i or nf (i+1) or not(gg (i-1) or gg i or gg(i+1) ) in
 119    for_all has_gflat (number mf);;
 120
 121 let flat_middle mf =
 122   let nf = nonflat_node mf in
 123   let has i = nf i or nf (i+1) or (part mf.edge i = E2 ) in
 124    for_all has (number mf);;
 125
 126 let minimal_node mf =
 127   let has i = (part mf.node i = N2) or (gminimal_edge mf i) or (gminimal_edge mf (i-1)) in
 128   for_all has (number mf);;
 129
 130 let minimal_node_flat mf =
 131   let has i = (nonflat_node mf i) or (part mf.node i = N2) or (gminimal_edge mf i && gminimal_edge mf (i-1)) in
 132   for_all has (number mf);;
 133
 134 let flat_extremal mf =
 135   let has i = nonflat_node mf i or nonflat_node mf (i+1) or bound_node mf i or bound_node mf (i+1) in
 136   for_all has (number mf);;
 137
 138 let extremal_node mf =
 139   let has i = flat_node mf i or flat_node mf (i+1) or flat_node mf (i+2) or bound_node mf (i+1) in
 140   for_all has (number mf);;
 141
 142 let flat_extremal_node mf =
 143    let has i = flat_node mf i or nonflat_node mf (i+1) or flat_node mf (i+2) or flat_node mf (i+3) or bound_node mf (i+1) or bound_node mf (i+2) in
 144   for_all has (number mf);;
 145
 146 let flat_extremal_node_sym mf =
 147    let has i = flat_node mf i or flat_node mf (i+1) or nonflat_node mf (i+2) or flat_node mf (i+3) or bound_node mf (i+1) or bound_node mf (i+2) in
 148   for_all has (number mf);;
 149
 150 let flat_count mf =
 151    length (filter not mf.nodeflats) > 2;;
 152
 153 let irreducible =
 154   fold_right (fun r s m -> r m && s m)
 155  [card;extreme_edge;flat_exists;no_triple_flat;balance;g_flat;
 156   flat_middle;minimal_node;minimal_node_flat;flat_extremal;
 157   extremal_node;flat_extremal_node;flat_extremal_node_sym;
 158   flat_count]
 159 (fun t->true);;
 160
 161
 162 (* symmetry reductions, add to the list of solutions only if it is shift-distinct from other solutions *)
 163
 164 let shift_one (x::xs) = xs @ [x];;
 165
 166 let shift a = { node = shift_one a.node; nodeflats = shift_one a.nodeflats; edge = shift_one a.edge };;
 167
 168 let rec shiftk k a =  if (k<=0) then a else shift (shiftk (k-1) a);;
 169
 170 let add_if_new a xs =
 171    if (exists (fun i -> mem (shiftk i a) xs) (upto (kn a))) then xs else a::xs;;
 172 let add_if_irred a xs = if (irreducible a ) then add_if_new a xs else xs;;
 173 (* let add_if_irred a xs = if (irreducible a ) then a ::xs else xs;; *)
 174
 175 let rec mk_all mkfn (imin,jmin,kmin) (imax,jmax,kmax) (i, j, k) acc  =
 176   let acc' = add_if_irred (mkfn (i, j, k)) acc in
 177   let mka = mk_all mkfn (imin,jmin,kmin) (imax,jmax,kmax) in
 178   if (k+1 < kmax) then mka (i,j,(k+1)) acc'
 179   else if (j+1 < jmax) then mka (i,(j+1),kmin) acc'
 180   else if (i+1 < imax) then mka ((i+1),jmin,kmin) acc'
 181   else acc';;
 182
 183 (* no flats *)
 184 let a3 = mk_all (mk_x 3 (3,1,3)) (0,0,0) (pow 3 3,1,pow 3 3) (0,0,0) [];;
 185 length a3;;
 186
 187 (* no flats *)
 188 let a4nf = mk_all (mk_x 4 (3,1,3)) (0,0,0) (pow 3 4,1,pow 3 4) (0,0,0) [];;
 189 length a4nf;;
 190
 191 (* exactly one flat *)
 192 let a4f = mk_all (mk_x 4 (3,2,3)) (0,1,0) (pow 3 4,2,pow 3 4) (0,1,0) [];;
 193 length a4f;;
 194
 195 (* no Gset *)
 196 (*
 197 let a5ng =mk_all (mk_x 5 (3,2,2)) (0,0,0) (pow 3 5,pow 2 5,pow 2 5) (0,0,0) [];;
 198 length a5ng;;
 199 *)
 200
 201 let a5 = mk_all (mk_x 5 (3,2,3)) (0,0,0) (pow 3 5,pow 2 5,pow 3 5) (0,0,0) [];;
 202 length a5;;
 203
 204 (* no Gset *)
 205 let a6 = mk_all (mk_x 6 (3,2,2)) (0,0,0) (pow 3 6,pow 2 6,pow 2 6) (0,0,0) [];;
 206 length a6;;
 207
 208 (* degrees of freedom *)
 209
 210 let freedom a = (kn a) - 3 + length (filter ((=) Nfree) a.node) - length (filter (fun t->t) a.nodeflats);;
 211 map freedom a4nf;;