Frame: optimize midpoint lookup

src/frame.c

@@ -117,7 +117,6 @@ int frame_lens_mesh_trace(FrameLensMesh *mesh, const SpacetimeSource *spacetime,
 typedef struct {
   size_t a, b, triangle;
   unsigned int side;
-  size_t midpoint;
 } MeshEdge;
 
 static int compare_mesh_edge(const void *left, const void *right) {
@@ -461,11 +460,18 @@ int frame_lens_mesh_finish_generation(FrameLensMesh *mesh,
   unsigned char *allowed = calloc(mesh->triangle_count, sizeof *allowed);
   signed char *parity = calloc(mesh->triangle_count, sizeof *parity);
   double *jacobians = calloc(mesh->triangle_count, sizeof *jacobians);
+  /* Sorted edge groups assign one midpoint to every incident triangle side.
+   * Keep that relation in the original triangle-side order so child emission
+   * stays O(T), rather than scanning every sorted edge for every child side. */
+  size_t *side_midpoints = malloc(edge_count * sizeof *side_midpoints);
   if (edges == NULL || requested == NULL || allowed == NULL || parity == NULL ||
-      jacobians == NULL) {
+      jacobians == NULL || side_midpoints == NULL) {
     free(edges); free(requested); free(allowed); free(parity); free(jacobians);
+    free(side_midpoints);
     return -1;
   }
+  for (size_t i = 0; i < edge_count; ++i)
+    side_midpoints[i] = SIZE_MAX;
   for (size_t i = 0; i < mesh->triangle_count; ++i) {
     const LensTriangle *triangle = &mesh->triangles[i];
     allowed[i] = triangle->level < config->max_level &&
@@ -473,7 +479,7 @@ int frame_lens_mesh_finish_generation(FrameLensMesh *mesh,
     for (unsigned int side = 0; side < 3; ++side) {
       size_t a = triangle->vertex[side], b = triangle->vertex[(side + 1) % 3];
       if (a > b) { size_t swap = a; a = b; b = swap; }
-      edges[3 * i + side] = (MeshEdge){a, b, i, side, SIZE_MAX};
+      edges[3 * i + side] = (MeshEdge){a, b, i, side};
     }
     if (allowed[i]) {
       const unsigned int side = longest_side(mesh, triangle);
@@ -568,6 +574,7 @@ int frame_lens_mesh_finish_generation(FrameLensMesh *mesh,
     mesh->sample_count = 0;
     mesh->samples_include_probes = 0;
     free(edges); free(requested); free(allowed); free(parity); free(jacobians);
+    free(side_midpoints);
     return 0;
   }
   /* Allocate a single stable midpoint vertex for each requested edge group. */
@@ -585,6 +592,7 @@ int frame_lens_mesh_finish_generation(FrameLensMesh *mesh,
   }
   if (ensure_vertices(mesh, mesh->vertex_count + midpoint_count)) {
     free(edges); free(requested); free(allowed); free(parity); free(jacobians);
+    free(side_midpoints);
     return -1;
   }
   size_t next_vertex = mesh->vertex_count;
@@ -604,7 +612,7 @@ int frame_lens_mesh_finish_generation(FrameLensMesh *mesh,
         midpoint = mesh->samples[probe_id].vertex;
       mesh->vertices[next_vertex++] = midpoint;
       for (size_t i = first; i < last; ++i)
-        edges[i].midpoint = next_vertex - 1;
+        side_midpoints[3 * edges[i].triangle + edges[i].side] = next_vertex - 1;
     }
     first = last;
   }
@@ -612,6 +620,7 @@ int frame_lens_mesh_finish_generation(FrameLensMesh *mesh,
   LensTriangle *children = calloc(old_count * 4, sizeof *children);
   if (children == NULL) {
     free(edges); free(requested); free(allowed); free(parity); free(jacobians);
+    free(side_midpoints);
     return -1;
   }
   size_t child_count = 0;
@@ -622,10 +631,7 @@ int frame_lens_mesh_finish_generation(FrameLensMesh *mesh,
     for (unsigned int side = 0; side < 3; ++side) {
       if (!requested[3 * t + side]) continue;
       ++count;
-      size_t a = triangle->vertex[side], b = triangle->vertex[(side + 1) % 3];
-      if (a > b) { size_t swap = a; a = b; b = swap; }
-      for (size_t e = 0; e < edge_count; ++e)
-        if (edges[e].a == a && edges[e].b == b) { middle[side] = edges[e].midpoint; break; }
+      middle[side] = side_midpoints[3 * t + side];
     }
     const size_t a = triangle->vertex[0], b = triangle->vertex[1], c = triangle->vertex[2];
     const unsigned int level = triangle->level + 1;
@@ -660,6 +666,7 @@ int frame_lens_mesh_finish_generation(FrameLensMesh *mesh,
   mesh->sample_count = 0;
   mesh->samples_include_probes = 0;
   free(edges); free(requested); free(allowed); free(parity); free(jacobians);
+  free(side_midpoints);
   return (int)midpoint_count;
 }