legacy/glpk/glpk/mkineq.ml

   1 (* create new inequalities for lp, cfsqp, formal spec *)
   2
   3 open Str;;
   4 open List;;
   5 open Num;;
   6 let sprintf = Printf.sprintf;;
   7
   8 (* from lpproc.ml *)
   9
  10 let unsplit d f = function
  11   | (x::xs) ->  fold_left (fun s t -> s^d^(f t)) (f x) xs
  12   | [] -> "";;
  13
  14 let join_lines  = unsplit "\n" (fun x-> x);;
  15
  16 let rec zip l1 l2 =
  17   match (l1,l2) with
  18         ([],[]) -> []
  19       | (h1::t1,h2::t2) -> (h1,h2)::(zip t1 t2)
  20       | _ -> failwith "zip";;
  21
  22
  23 (* end from lpproc.ml *)
  24
  25 type sgn = GT | GE;;
  26
  27 type constrain = NONE | SMALLTRI | BIGTRI ;;
  28
  29 type node = LOW | MEDIUM | HIGH | EXTRA ;;
  30
  31 type edge = SMALL | BIG;;
  32
  33 type decimal = Dec of string | Sqrt2 | Sqrt8;;
  34
  35 let has_point s = try
  36   let _ = String.index s '.' in true
  37   with failure -> false;;
  38
  39 let add_point s = if has_point s then s else s^".0";;
  40
  41 let dec_of_string s =
  42   if s="s2" then Sqrt2
  43   else if s="s8" then Sqrt8
  44   else Dec s;;
  45
  46 let d = dec_of_string;;
  47
  48 let ds:string -> decimal list = fun s ->
  49   let ss = split_sp s in
  50     map dec_of_string ss;;
  51
  52 let float_of_dec:decimal -> float = function
  53   | Sqrt2 -> sqrt(2.0)
  54   | Sqrt8 -> sqrt(8.0)
  55   | Dec x -> float_of_string x;;
  56
  57 let holtext_of_dec:decimal -> string = function
  58   | Sqrt2 -> "s2"
  59   | Sqrt8 -> "s8"
  60   | Dec x ->  let y = add_point x in
  61       if (y.[0]= '-') then " -- #"^(String.sub y 1 (String.length y - 1))
  62       else "#"^y;;
  63
  64
  65 holtext_of_dec (Dec "2.0");;
  66
  67 let holtext_of_declist:decimal list -> string = fun xs ->
  68   "["^(unsplit ";" holtext_of_dec xs)^"]";;
  69
  70 let holtext_of_sgn:sgn -> string = function
  71   | GT -> "( > )"
  72   | GE -> "( >= )";;
  73
  74 let holtext_of_constrain:constrain -> string = function
  75   | NONE -> "Cnone"
  76   | BIGTRI -> "Cstd3_big"
  77   | SMALLTRI -> "Cstd3_small";;
  78
  79 let split_sp=  Str.split (regexp " +");;
  80
  81
  82
  83 (*
  84 represents on rectangle xmin[6], xmax[6]:
  85
  86 azim[i]*azim i y + rhzim[i]*rhzim i y + tau0 * taumar y + sol0 * sol y
  87 sgn
  88 c0 + c dot (y-p).
  89
  90 *)
  91
  92 type ineq = {
  93   mutable id : string;
  94   mutable constrain: constrain;
  95   mutable sgn : sgn;
  96   mutable xmin : decimal list;
  97   mutable xmax : decimal list;
  98   mutable c0 : decimal;
  99   mutable c : decimal list;
 100   mutable p : decimal list;
 101   mutable azim : decimal list;
 102   mutable rhzim : decimal list;
 103   mutable tau0 : decimal;
 104   mutable sol0 : decimal;
 105 };;
 106
 107 let hh =  {
 108   id = "21444";
 109   constrain = NONE;
 110   sgn = GT;
 111   xmin = ds "2 2 2 2 2 2";
 112   xmax = ds "2.52 2.52 2.52 2.52 2.52 2.52";
 113   c0 = d "2.34";
 114   c = ds "2.0 2.1 2.2 2.3 2.4 2.5";
 115   p = ds "3.0 3.1 3.2 3.3 3.4 3.5";
 116   azim = ds "1 1 1";
 117   rhzim = ds "-1.0 -2.0 -3.0";
 118   tau0 = d "1.0";
 119   sol0 = d "-4.0";
 120 };;
 121
 122
 123 let holtext_of_ineq:ineq->string = fun h ->
 124    let p = sprintf in
 125      join_lines [
 126        p"let hol_ineq%s = `hol_ineq (\"%s\", " h.id h.id;
 127        p"  %s," (holtext_of_constrain h.constrain);
 128        p"  %s," (holtext_of_sgn h.sgn);
 129        p"  %s," (holtext_of_declist h.xmin);
 130        p"  %s," (holtext_of_declist h.xmax);
 131        p"  %s," (holtext_of_dec h.c0);
 132        p"  %s," (holtext_of_declist h.c);
 133        p"  %s," (holtext_of_declist h.p);
 134        p"  %s," (holtext_of_declist h.azim);
 135        p"  %s," (holtext_of_declist h.rhzim);
 136        p"  %s," (holtext_of_dec h.tau0);
 137        p"  %s)` " (holtext_of_dec h.sol0);
 138      ];;
 139 holtext_of_ineq hh;;
 140 (* ampl text generation of triangle ineqs
 141    ocaml numbering 012345
 142    ampl numbering 123456 *)
 143
 144
 145 let aug i = i+1;;
 146 let nz s = (float_of_dec s <> 0.0);;
 147 let hasnz s = exists nz s;;
 148 let unempty   = filter (fun t -> t <> "");;
 149
 150 let ampl_of_dec:decimal -> string = function
 151   | Sqrt2 -> "+1.4142135623730951"
 152   | Sqrt8 -> "+2.8284271247461903"
 153   | Dec x ->
 154       if (x.[0]= '-') then x else "+"^x;;
 155
 156 let comp:decimal->decimal->int =
 157   fun a b ->
 158     if (a=b) then 0 else compare (float_of_dec a) (float_of_dec b);;
 159
 160 let less_equal bs cs =
 161   []= filter (fun t -> comp (fst t) (snd t) >0 ) (zip bs cs);;
 162
 163 let domain_covers (lo,high) h =
 164   less_equal h.xmin lo && less_equal high h.xmax;;
 165
 166 let domain_covers_itriangle =
 167   domain_covers (ds "2 2 2 2 2 2",ds "2.52 2.52 2.52 2.52 2.52 2.52");;
 168
 169 let domain_covers_apiece =
 170   domain_covers (ds "2 2 2 2 2.52 2.52",ds "2.52 2.52 2.52 2.52 s8 s8");;
 171
 172 let domain_covers_flat =
 173   domain_covers (ds "2 2 2 2.52 2 2",ds "2.52 2.52 2.52 s8 2.52 2.52");;
 174
 175 let domain_covers_apex_sup_flat =
 176   domain_covers (ds "2 2 2 s8 2 2",ds "2.52 2.52 2.52 3.0 2.52 2.52");;
 177
 178 let domain_covers_std3_big h =
 179   domain_covers_itriangle h && (h.constrain = BIGTRI);;
 180
 181 let domain_covers_std3_small h =
 182   domain_covers_itriangle h && (h.constrain = SMALLTRI);;
 183
 184 let node_range = function
 185   | LOW -> ds "2 2.18"
 186   | MEDIUM -> ds "2.18 2.36"
 187   | HIGH -> ds "2.18 2.52"
 188   | EXTRA -> ds "2.36 2.52";;
 189
 190 let domain_covers_f f node h i =
 191   let ymin = nth h.xmin i in
 192   let ymax = nth h.xmax i in
 193   let [mm;mx] = f node in
 194     less_equal [ymin] [mm] && less_equal [mx] [ymax];;
 195
 196 let domain_covers_node = domain_covers_f node_range;;
 197
 198 let edge_range = function
 199   | SMALL -> ds "2 2.25"
 200   | BIG -> ds "2.25 2.52";;
 201
 202 let domain_covers_edge = domain_covers_f edge_range;;
 203
 204 let string_of_domain h = "";;
 205
 206 let ampltext_of_ineq:ineq->string = fun h ->
 207   let p = sprintf in
 208   let a = ampl_of_dec in
 209   let mkone f s = if nz f then p"  %s * %s " (a f) s else "" in
 210   let mk_madd j = p"  %s * (y%d[i2,j] - (%s))"
 211     (a (nth h.c j)) (j+1) (a (nth h.p j)) in
 212     join_lines (unempty[
 213       "# ";
 214       p"ineq%s 'ID[%s]' {(i1,i2,i3,j) in e_dart : " h.id h.id;
 215       string_of_domain h;
 216       "}:";
 217       mkone h.tau0 "tau[j]";
 218       mkone h.sol0 "sol[j]";
 219       mkone (nth h.azim 0) "azim[i1,j]";
 220       mkone (nth h.azim 1) "azim[i2,j]";
 221       mkone (nth h.azim 2) "azim[i3,j]";
 222       mkone (nth h.rhzim 0) "rhzim[i1,j]";
 223       mkone (nth h.rhzim 1) "rhzim[i2,j]";
 224       mkone (nth h.rhzim 2) "rhzim[i3,j]";
 225       " >= ";
 226       p"  %s" (a h.c0);
 227       mk_madd 0;mk_madd 1;mk_madd 2;mk_madd 3;mk_madd 4;mk_madd 5
 228     ]);;
 229 1;;
 230 ampltext_of_ineq hh;;