6 Notes on the paper "Contact Graphs of Unit Sphere Packings Revisited"
8 preprint date October 17, 2012.
11 We show the following results:
13 There are no kissing arrangements of 12 spheres with at least 25 contact edges.
14 There are no kissing arrangements of 12 spheres with at least 11 contact triples.
17 This is not a stand-alone piece of code.
18 For the graph generation, it requires the installation of java.
20 It also requires parts of the Flyspeck project, especially
21 Graph_control for graph generation, and Glpk_link to process the raw data,
23 We also use some functions such as (--) and chop_list from HOL-Light's lib.ml.
25 As the code is currently written, it assumes that the entire Flyspeck project
26 has been loaded, but with a bit of work, it could be made to depend on just
27 the indicated modules.
31 The relevant output is coord_edge and coord_triple, which exhibit inconsistent
32 coordinate assignments for the nodes of the graphs.
36 flyspeck_needs ("../graph_generator/graph_control.hl");;
37 flyspeck_needs "../glpk/glpk_link.ml";;
40 module Bezdek_reid = struct
46 (* The file names need to be adjusted for different users *)
47 let modelfile = "/Users/thomashales/Desktop/googlecode/flyspeck/projects_discrete_geom/bezdek_reid/br_model.mod";;
49 (* This is the default output location for the graph generation, as set in the java code. *)
50 let tmpfile = "/tmp/graph_out.txt";;
52 let bezdek_reid_properties_edge =
55 "K. Bezdek and S. Reid: Contact Graphs of Unit Sphere Packings Revisited (2012)";
57 exclude_degree2=false;
58 exclude_pent_qrtet=true;
59 exclude_2_in_quad=true;
60 exclude_1_in_tri=true;
62 (* require exactly 12 spheres *)
66 (* degree at most 5 at every node *)
68 node_card_max_at_exceptional_vertex=5;
70 (* if the are more than 5.5 edges removed, then ignore *)
74 (* a quad face = 1 edge removed, pent = 2 edges removed, etc. *)
75 table_weight_d = [0;0;0;0;10;20;30;40;50];
76 table_weight_a = [(0,0);(1,0);(2,0);(3,0);(4,0)];
77 table_weight_b = [(0,3,30);(0,4,40);(0,5,50);
84 let bezdek_reid_properties_triplet =
87 "K. Bezdek and S. Reid: Contact Graphs of Unit Sphere Packings Revisited (2012)";
89 exclude_degree2=false;
90 exclude_pent_qrtet=true;
91 exclude_2_in_quad=true;
92 exclude_1_in_tri=true;
94 (* require exactly 12 spheres *)
98 (* degree at most 5 at every node *)
100 node_card_max_at_exceptional_vertex=5;
102 (* if the are more than 9.5 non-contact edges, then ignore *)
106 (* a quad face = 1 edge removed, pent = 2 edges removed, etc. *)
107 (* allow up to an 11-gon. 11-gon, 10-gon, 9-gon ruled out
108 by extra n-gons at remaining vertices. *)
109 table_weight_d = [0;0;0;0;20;30;40;50;60];
110 table_weight_a = [(0,0);(1,0);(2,0);(3,0);(4,0)];
111 table_weight_b = [(0,3,60);(0,4,80);(0,5,100);
119 let max_degree hypl =
120 let fl = List.flatten (map (fun h-> ((mk_bb h).std_faces)) hypl) in
121 let ll = map List.length fl in
124 (* initial triangle *)
125 let sqrt = Pervasives.sqrt;;
126 let p1 = (2.0,0.0,0.0);;
127 let p2 = (1.0,sqrt(3.0),0.0);;
128 let p3 = (1.0,1.0/. sqrt(3.0), sqrt(8.0) /. sqrt(3.0));;
131 (* We include a few vector functions from tikz.ml *)
133 let (+...) (x1,x2,x3) (y1,y2,y3) = (x1 +. y1, x2 +. y2, x3+. y3);;
135 let (-...) (x1,x2,x3) (y1,y2,y3) = (x1 -. y1, x2 -. y2, x3-. y3);;
137 let ( %... ) s (x1,x2,x3) = (s *. x1, s *. x2, s*. x3);;
139 let ( *... ) (x1,x2,x3) (y1,y2,y3) = (x1 *. y1 +. x2 *. y2 +. x3 *. y3);;
141 let cross (x1,x2,x3) (y1,y2,y3) =
142 (x2 *. y3 -. x3 *. y2, x3 *. y1 -. x1 *. y3, x1 *. y2 -. x2 *. y1);;
144 let normalize3 x = (1.0 /. sqrt(x *... x)) %... x;;
147 let z = x -... y in sqrt (z *... z);;
150 let u = normalize3 (cross b c) in
151 a -... ((2.0 *. (u *... a)) %... u);;
153 let list_of (a,b,c) = [a;b;c];;
155 let tuple_of [a;b;c] = (a,b,c);;
157 let common a b = intersect (list_of a) (list_of b);;
159 let adj3 a b = (length (list_of a) = 3) && (length (list_of b) = 3) &&
160 length (common a b) = 2;;
165 | a::r -> (map (fun i -> (a,i)) y) @ (outer r y);;
169 [] -> failwith "find"
170 | (h::t) -> if p(h) then h else find p t;;
172 let find_adj_pair al bl =
173 let (s,t) = find (fun (s,t) -> adj3 s t) (outer al bl) in
174 let c = common s t in
175 let (c1,c2) = (List.nth c 0,List.nth c 1) in
176 let (sl,tl) = (list_of s, list_of t) in
177 let a = hd (subtract sl c) in
178 let b = hd (subtract tl c) in
181 let update_coord (cl,br1 ,br2) =
182 let (a,b,c,a',s,t) = find_adj_pair br1 br2 in
183 let pa = assoc a cl in
184 let pb = assoc b cl in
185 let pc = assoc c cl in
186 let pa' = reflect pa pb pc in
187 ((a',pa')::cl , t::br1 , subtract br2 [t]);;
189 let rec mk_coords (cl,br1,br2) =
190 if (br2 = [] or not( can (update_coord) (cl,br1,br2)))
191 then (cl,br1,br2) else mk_coords (update_coord (cl,br1,br2));;
193 let triangles_in_archive_string ll i =
194 let case1 = List.nth ll i in
195 let br = (mk_bb (case1)).std_faces in
196 let br_tri = filter (fun t -> List.length t = 3) br in
197 map tuple_of br_tri;;
199 let do_one_case brt =
200 let (br1,br2) = chop_list 1 (brt) in
201 let [(a1,a2,a3)] = br1 in
202 let coordlist = [(a1,p1);(a2,p2);(a3,p3)] in
203 mk_coords (coordlist,br1,br2);;
205 let rec mk_all_coords (cll,br2) =
206 if br2 = [] then (cll,[])
208 let (cl,_,br2') = do_one_case br2 in
209 mk_all_coords (cl::cll, br2');;
211 (* consistency checks on coordinates *)
214 let vl = setify (map fst dc) in
215 let dis_ok = map (fun (i,j) -> if (j<=i) then true
216 else (1.99 < dist3 (assoc i dc) (assoc j dc)))
218 (setify dis_ok = [true]);;
220 let self_consis_one dc i =
222 let dci = filter (fun (j,_) -> (j =i)) dc in
223 let dr = map snd dci in
224 let dis_ok = map (fun (a,b) -> dist3 a b < eps) (outer dr dr) in
225 setify dis_ok = [true];;
228 let vl = setify (map fst dc) in
229 forall (self_consis_one dc) vl;;
231 let inter_consis_one dc1 dc2 =
232 let vl = intersect (map fst dc1) (map fst dc2) in
233 let a1 i = assoc i dc1 in
234 let a2 i = assoc i dc2 in
235 let abs = Pervasives.abs_float in
237 forall (fun (i,j) -> abs(dist3 (a1 i) (a1 j) -. dist3 (a2 i) (a2 j)) < eps)
240 let inter_consis cll =
241 forall (fun (a,b) -> inter_consis_one a b) (outer cll cll);;
244 let br = triangles_in_archive_string ll n in
245 let (cll,_) = mk_all_coords ([],br) in
248 let test_coords ll n =
249 let cll = mk_c ll n in
250 (forall (is_packing) cll) &&
251 (forall (self_consis) cll) &&
254 (* we see by inspection from the explicit coordinates that they are
255 overdetermined. test_edge and test_triplet should both evaluate to [false] *)
257 (* edge case processing *)
260 let archiveraw_edge = tmpfile in
261 let _ = Graph_control.run bezdek_reid_properties_edge in
262 Glpk_link.strip_archive archiveraw_edge;;
264 let count_edge = List.length hypermap_edge;; (* 18 cases *)
265 let max_degree_edge = max_degree hypermap_edge;;
266 let range_edge = (0--(List.length hypermap_edge - 1));;
267 let test_edge = setify (map (test_coords hypermap_edge) range_edge);;
272 let hypermap_triplet =
273 let archiveraw_triplet = tmpfile in
274 let _ = Graph_control.run bezdek_reid_properties_triplet in
275 Glpk_link.strip_archive archiveraw_triplet;;
276 let count_triplet = List.length hypermap_triplet;; (* 1335 cases *)
277 let max_degree_edge = max_degree hypermap_triplet;; (* 8 *)
278 let range_triplet = (0--(List.length hypermap_triplet - 1));;
279 let test_triplet = setify (map (test_coords hypermap_triplet) range_triplet);;