From 13f94397661f67256a68849da13d2029880647c2 Mon Sep 17 00:00:00 2001
From: Ian Jackson <ian@davenant.greenend.org.uk>
Date: Sun, 6 Jan 2008 20:59:54 +0000
Subject: [PATCH] quite nice

---
 energy.c   | 51 +++++++++++++++++++++++++++++++++++++++------------
 minimise.h |  1 +
 2 files changed, 40 insertions(+), 12 deletions(-)

diff --git a/energy.c b/energy.c
index 49d4e43..4723e0c 100644
--- a/energy.c
+++ b/energy.c
@@ -27,9 +27,10 @@ double compute_energy(const struct Vertices *vs) {
 
   if (printing) printf("cost > energy |");
 
-  COST(1e2, edgewise_vertex_displacement_cost(vs->a));
-  COST(1e2, graph_layout_cost(vs->a));
+  COST(3e2, edgewise_vertex_displacement_cost(vs->a));
   COST(1e3, edge_length_variation_cost(vs->a));
+  COST(0.2e3, rim_proximity_cost(vs->a));
+//  COST(1e2, graph_layout_cost(vs->a));
   COST(1e8, noncircular_rim_cost(vs->a));
 
   if (printing) printf("| total %# e |", energy);
@@ -57,7 +58,7 @@ double compute_energy(const struct Vertices *vs) {
 
 static void addcost(double *energy, double tweight, double tcost, int pr) {
   double tenergy= tweight * tcost;
-  if (pr) printf(" %# e > %# e |", tcost, tenergy);
+  if (pr) printf(" %# e > %# e* |", tcost, tenergy);
   *energy += tenergy;
 }
 
@@ -181,22 +182,48 @@ double edge_length_variation_cost(const Vertices vertices) {
   return cost;
 }    
 
+/*---------- rim proximity cost ----------*/
+
+static void find_nearest_oncircle(double oncircle[D3], const double p[D3]) {
+  /* By symmetry, nearest point on circle is the one with
+   * the same angle subtended at the z axis. */
+  oncircle[0]= p[0];
+  oncircle[1]= p[1];
+  oncircle[2]= 0;
+  double mult= 1.0/ magnD(oncircle);
+  oncircle[0] *= mult;
+  oncircle[1] *= mult;
+}
+
+double rim_proximity_cost(const Vertices vertices) {
+  double oncircle[3], cost=0;
+  int v;
+
+  FOR_VERTEX(v) {
+    int y= v >> YSHIFT;
+    int nominal_edge_distance= y <= Y/2 ? y : Y-1-y;
+    if (nominal_edge_distance==0) continue;
+
+    find_nearest_oncircle(oncircle, vertices[v]);
+
+    cost +=
+      vertex_mean_edge_lengths[v] *
+      (nominal_edge_distance*nominal_edge_distance) /
+      (hypotD2(vertices[v], oncircle) + 1e-6);
+  }
+  return cost;
+}
+
 /*---------- noncircular rim cost ----------*/
 
 double noncircular_rim_cost(const Vertices vertices) {
   int vy,vx,v;
   double cost= 0.0;
+  double oncircle[3];
 
   FOR_RIM_VERTEX(vy,vx,v) {
-    double oncircle[3];
-    /* By symmetry, nearest point on circle is the one with
-     * the same angle subtended at the z axis. */
-    oncircle[0]= vertices[v][0];
-    oncircle[1]= vertices[v][1];
-    oncircle[2]= 0;
-    double mult= 1.0/ magnD(oncircle);
-    oncircle[0] *= mult;
-    oncircle[1] *= mult;
+    find_nearest_oncircle(oncircle, vertices[v]);
+    
     double d2= hypotD2(vertices[v], oncircle);
     cost += d2*d2;
   }
diff --git a/minimise.h b/minimise.h
index 2043d16..f1da76f 100644
--- a/minimise.h
+++ b/minimise.h
@@ -19,6 +19,7 @@ extern double vertex_areas[N], vertex_mean_edge_lengths[N], edge_lengths[N][V6];
 double noncircular_rim_cost(const Vertices vertices);
 double edgewise_vertex_displacement_cost(const Vertices vertices);
 double edge_length_variation_cost(const Vertices vertices);
+double rim_proximity_cost(const Vertices vertices);
 
 extern const char *input_file, *output_file;
 extern char *output_file_tmp;
-- 
2.30.2