math-research/papers/level_switching/experiments/leaf_distance_algorithm.py

"""User-proposed simpler algorithm:

1) Assign each face x = min dual-tree distance to any leaf face (ear).
2) When face F at depth d has no balanced surface switch, edge-switch
   on the edge between F and F's neighbour with the smallest x.
3) Repeat until the number of depth-d faces has strictly decreased.
"""
import os
import sys
sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)))
import math
import networkx as nx
from stress_test_termination import (
    compute_depths, apply_switch, check_balanced, face_edges
)


def compute_x(faces, outer_set):
    """x(F) = min dual-tree distance from F to any leaf face."""
    D = nx.Graph()
    D.add_nodes_from(range(len(faces)))
    for i, fi in enumerate(faces):
        for j, fj in enumerate(faces):
            if i < j and face_edges(fi) & face_edges(fj):
                D.add_edge(i, j)
    leaves = [i for i in D.nodes if D.degree(i) == 1]
    x = {}
    for i in range(len(faces)):
        if not leaves:
            x[i] = float('inf')
            continue
        x[i] = min(nx.shortest_path_length(D, i, lf) for lf in leaves)
    return x, D


def neighbour_at_edge(faces, F_idx, e):
    for j, fj in enumerate(faces):
        if j != F_idx and e in face_edges(fj):
            return j
    return None


def run_leaf_distance_algorithm(faces, outer_set, max_steps=500,
                                verbose=True):
    """Run the algorithm. If F admits a balanced switch, take it; else
    pick the neighbour with smallest x and switch on that edge."""
    total = 0
    log = []
    for step in range(max_steps):
        depth = compute_depths(faces, outer_set)
        d_max = max(depth.values())
        if d_max == 0:
            log.append(f'terminated at step {step}')
            return faces, total, log

        max_d_faces = [i for i, d in depth.items() if d == d_max]
        F_idx = max_d_faces[0]
        F = faces[F_idx]

        ok, _, fp_idx, e = check_balanced(F_idx, faces, depth, outer_set)
        if ok:
            Fp = faces[fp_idx]
            u, v = tuple(e)
            w = [vert for vert in F if vert not in (u, v)][0]
            x = [vert for vert in Fp if vert not in (u, v)][0]
            log.append(f'step {step}: BALANCED on edge ({u},{v}) '
                       f'with F\' = {Fp}')
            faces = apply_switch(faces, (u, v), (w, x))
            total += 1
            continue

        # Find x values
        x_vals, D = compute_x(faces, outer_set)
        # Among neighbours of F (interior edges, hence in dual graph),
        # pick the one with smallest x.
        nb_choices = []
        for e_test in face_edges(F):
            if e_test in outer_set:
                continue
            nb_idx = neighbour_at_edge(faces, F_idx, e_test)
            if nb_idx is None:
                continue
            nb_choices.append((x_vals[nb_idx], e_test, nb_idx))
        if not nb_choices:
            log.append(f'step {step}: no interior neighbour; stuck')
            break
        nb_choices.sort()
        x_min, e_chosen, nb_idx = nb_choices[0]
        Fnb = faces[nb_idx]
        u, v = tuple(e_chosen)
        w = [vert for vert in F if vert not in (u, v)][0]
        xv = [vert for vert in Fnb if vert not in (u, v)][0]
        log.append(f'step {step}: x-preprocess on edge ({u},{v}) '
                   f'(F\'={Fnb}, x={x_min}; other choices: '
                   f'{[(c[0]) for c in nb_choices[1:]]})')
        faces = apply_switch(faces, (u, v), (w, xv))
        total += 1

    log.append(f'hit max_steps; final max depth = '
               f'{max(compute_depths(faces, outer_set).values())}')
    return faces, total, log


# ----- TEST CASE 1: 9-vertex example -----
n9 = 9
outer9 = [(i, (i + 1) % n9) for i in range(n9)]
outer_set9 = {frozenset(e) for e in outer9}
faces9 = [
    (0, 1, 2), (0, 2, 3), (3, 4, 5), (3, 5, 6),
    (6, 7, 8), (6, 8, 0), (0, 3, 6),
]
print('=== 9-vertex example ===')
_, total, log = run_leaf_distance_algorithm(faces9, outer_set9)
print('\n'.join(log))
print(f'total switches: {total}\n')


# ----- TEST CASE 2: 24-vertex doubly-lopsided example -----
n24 = 24
outer24 = [(i, (i + 1) % n24) for i in range(n24)]
outer_set24 = {frozenset(e) for e in outer24}

def arm(a, b):
    return [
        (a, a + 1, a + 2), (a, a + 2, b), (a + 2, a + 3, a + 4),
        (a + 2, a + 4, b), (a + 4, a + 5, a + 6), (a + 4, a + 6, b),
        (a + 6, a + 7, b),
    ]
faces24 = [(0, 8, 16)]
for (a, b) in [(0, 8), (8, 16), (16, 24)]:
    fs = arm(a, b)
    fs = [tuple(0 if v == 24 else v for v in vt) for vt in fs]
    faces24.extend(fs)
print('=== 24-vertex doubly-lopsided example ===')
_, total, log = run_leaf_distance_algorithm(faces24, outer_set24)
print('\n'.join(log[:25]))
if len(log) > 25:
    print(f'... ({len(log) - 25} more lines)')
print(f'total switches: {total}\n')


# ----- TEST CASE 3: random outerplanar n=40 (one of the previously-slow seeds) -----
from stress_test_termination import random_triangulation
seed = 94476710
outer, _, faces = random_triangulation(40, seed=seed)
outer_set = {frozenset(e) for e in outer}
print(f'=== Random n=40 (seed={seed}) with the new algorithm ===')
_, total, log = run_leaf_distance_algorithm(faces, outer_set, max_steps=1000)
print(f'total switches: {total} (random algorithm with cap=10000 needed 863)')
print(f'first 5 lines: {log[:5]}')
print(f'last 2 lines: {log[-2:]}')