X-Git-Url: http://www.chiark.greenend.org.uk/ucgi/~ian/git?a=blobdiff_plain;f=bgl.cpp;h=193a75a97fac5eac1b901aa2afcb7ff511a572cf;hb=8dd86f878d8fd02ebcda87454065e3bb18d13432;hp=4b06897df09a4ba5e65a903427ab021ce9261cc9;hpb=83051837693783280938bf36e9e43c7687736753;p=moebius2.git

diff --git a/bgl.cpp b/bgl.cpp
index 4b06897..193a75a 100644
--- a/bgl.cpp
+++ b/bgl.cpp
@@ -22,10 +22,11 @@ extern "C" {
 }
 
 /*
- * edge descriptor f =   00 | e | y     | x
- *                            3  YBITS   XBITS
+ * edge descriptor f =   0000 | e | y     | x
+ *                              3  YBITS   XBITS
  *
- * e is 0..5.  The edge is edge e out of vertex (x,y).
+ * e is 0..6.  The edge is edge e out of vertex (x,y), or if
+ * e==6 it's the `at end' value for the out edge iterator.
  *
  * BGL expects an undirected graph's edges to have two descriptors
  * each, one in each direction (otherwise e would be just 0..2).
@@ -46,8 +47,6 @@ extern "C" {
 #define VMASK (YMASK|XMASK)
 #define ESHIFT (YBITS+XBITS)
 
-class Graph { }; // this is a dummy as our graph has no actual representation
-
 using namespace boost;
 
 /*
@@ -93,17 +92,19 @@ class OutEdgeIterator :
   OutEdgeIterator() { }
   OutEdgeIterator(int _f) : f(_f) { }
   OutEdgeIterator(int v, int e) : f(e<<ESHIFT | v) {
-    //printf("constructed v=%x e=%x f=%03x\n",v,e,f);
+    //printf("constructed v=%02x e=%x f=%03x\n",v,e,f);
   }
 
   static int voe_min(int _v) { return (_v & YMASK) ? 2 : 3; }
-  static int voe_max(int _v) { return (~_v & YMASK) ? V6 : 4; }
-  static int voe_degree(int _v) { return (_v & YMASK | ~_v & YMASK) ? 4 : V6; }
+  static int voe_max(int _v) { return (_v & YMASK)==(Y-1) ? V6 : 4; }
+  static int voe_degree(int _v) { return RIM_VERTEX_P(_v) ? 4 : V6; }
 };
  
 typedef counting_iterator<int> VertexIterator;
 
 namespace boost {
+  class Graph { }; // this is a dummy as our graph has no actual representation
+
   // We make Graph a model of various BGL Graph concepts.
   // This mainly means that graph_traits<Graph> has lots of stuff.
 
@@ -114,6 +115,7 @@ namespace boost {
     public virtual vertex_list_graph_tag,
     public virtual edge_list_graph_tag { };
 
+  template<>
   struct graph_traits<Graph> {
     // Concept Graph:
     typedef int vertex_descriptor; /* vertex number, -1 => none */
@@ -136,7 +138,11 @@ namespace boost {
 
   // Concept IncidenceGraph:
   inline int source(int f, const Graph&) { return f&VMASK; }
-  inline int target(int f, const Graph&) { return EDGE_END2(f&VMASK, f>>ESHIFT); }
+  inline int target(int f, const Graph&) {
+    int v2= EDGE_END2(f&VMASK, f>>ESHIFT);
+    //printf("traversed %03x..%02x\n",f,v2);
+    return v2;
+  }
   inline std::pair<OutEdgeIterator,OutEdgeIterator>
   out_edges(int v, const Graph&) {
     return std::make_pair(OutEdgeIterator(v, OutEdgeIterator::voe_min(v)),
@@ -147,7 +153,8 @@ namespace boost {
   }
 
   // Concept VertexListGraph:
-  inline std::pair<VertexIterator,VertexIterator> vertices(const Graph&) {
+  inline
+  std::pair<VertexIterator,VertexIterator> vertices(const Graph&) {
     return std::make_pair(VertexIterator(0), VertexIterator(N));
   }
   inline unsigned num_vertices(const Graph&) { return N; }
@@ -183,7 +190,7 @@ double graph_layout_cost(const Vertices v, const double vertex_areas[N]) {
    */
   static const double d2_epsilon= 1e-6;
   
-  double edge_weights[N*V6], vertex_distances[N], total_cost=0;
+  double edge_weights[V6<<ESHIFT], vertex_distances[N], total_cost=0;
   int v1,v2,e,f;
 
   FOR_VERTEX(v1)
@@ -196,11 +203,19 @@ double graph_layout_cost(const Vertices v, const double vertex_areas[N]) {
     double a1= vertex_areas[v1];
     single_source_shortest_paths(v1, edge_weights, vertex_distances);
     FOR_VERTEX(v2) {
+      if (v1 == v2) continue;
       double a2= vertex_areas[v2];
       double d2= hypotD2plus(v[v1],v[v2], d2_epsilon);
-      double sd= vertex_distances[v2] / d2;
-      double sd2= sd*sd;
-      total_cost += a1*a2 * (sd2 - 1) / (d2*d2);
+      double s= vertex_distances[v2];
+      double s2= s*s + d2_epsilon;
+      double sd2= s2 / d2;
+      double cost_contrib= a1*a2 * (sd2 - 1) / (d2*d2);
+      if (cost_contrib < -1e-4) {
+	printf("layout %03x..%03x (a=%g,%g) s=%g s2=%g d2=%g sd2=%g"
+	       " cost+=%g\n", v1,v2, a1,a2, s,s2,d2,sd2, cost_contrib);
+	abort();
+      }
+      total_cost += cost_contrib;
     }
   }
   return total_cost;