diff --git a/.gitattributes b/.gitattributes new file mode 100644 index 0000000..5c0d01e --- /dev/null +++ b/.gitattributes @@ -0,0 +1 @@ +docs/unweighted-sat-search-report.pdf binary diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml index ea5893d..059a13a 100644 --- a/.github/workflows/CI.yml +++ b/.github/workflows/CI.yml @@ -28,7 +28,6 @@ jobs: os: - ubuntu-latest - macos-latest - - windows-latest arch: - x64 steps: @@ -76,4 +75,4 @@ jobs: using Documenter: DocMeta, doctest using GadgetSearch DocMeta.setdocmeta!(GadgetSearch, :DocTestSetup, :(using GadgetSearch); recursive=true) - doctest(GadgetSearch) \ No newline at end of file + doctest(GadgetSearch) diff --git a/Project.toml b/Project.toml index ae5b2e3..857636f 100644 --- a/Project.toml +++ b/Project.toml @@ -16,6 +16,7 @@ IterTools = "c8e1da08-722c-5040-9ed9-7db0dc04731e" JSON3 = "0f8b85d8-7281-11e9-16c2-39a750bddbf1" JuMP = "4076af6c-e467-56ae-b986-b466b2749572" Karnak = "cd156443-31ad-4f6f-850f-a93ee5f75905" +Kissat_jll = "971a2a11-349f-5eb6-a4bf-bc548fa92b3a" Luxor = "ae8d54c2-7ccd-5906-9d76-62fc9837b5bc" ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" @@ -34,6 +35,7 @@ IterTools = "1.10.0" JSON3 = "1.14.1" JuMP = "1" Karnak = "1.1.0" +Kissat_jll = "3.0.0" Luxor = "4.1.0" ProgressMeter = "1.10.4" Serialization = "1.11.0" diff --git a/README.md b/README.md index a693bf5..68a9878 100644 --- a/README.md +++ b/README.md @@ -47,6 +47,8 @@ Unweighted mode ignores vertex weights, treats the pins as open vertices, and computes the reduced alpha tensor. The lattice dimensions can be adjusted from 2×2 to 20×20 and are preserved in exported JSON. +The unweighted SAT search uses Kissat and supports Linux and macOS. + ### Rydberg Model (MIS-based) ```julia diff --git a/docs/figures/unweighted-sat/bdd-example.svg b/docs/figures/unweighted-sat/bdd-example.svg new file mode 100644 index 0000000..91f2cfd --- /dev/null +++ b/docs/figures/unweighted-sat/bdd-example.svg @@ -0,0 +1 @@ + \ No newline at end of file diff --git a/docs/figures/unweighted-sat/connectivity-layers.svg b/docs/figures/unweighted-sat/connectivity-layers.svg new file mode 100644 index 0000000..897841a --- /dev/null +++ b/docs/figures/unweighted-sat/connectivity-layers.svg @@ -0,0 +1 @@ + \ No newline at end of file diff --git a/docs/figures/unweighted-sat/frame-enumeration.svg b/docs/figures/unweighted-sat/frame-enumeration.svg new file mode 100644 index 0000000..ba48cb3 --- /dev/null +++ b/docs/figures/unweighted-sat/frame-enumeration.svg @@ -0,0 +1 @@ + \ No newline at end of file diff --git a/docs/figures/unweighted-sat/pipeline.svg b/docs/figures/unweighted-sat/pipeline.svg new file mode 100644 index 0000000..3b0566c --- /dev/null +++ b/docs/figures/unweighted-sat/pipeline.svg @@ -0,0 +1 @@ + \ No newline at end of file diff --git a/docs/figures/unweighted-sat/verification.svg b/docs/figures/unweighted-sat/verification.svg new file mode 100644 index 0000000..d272968 --- /dev/null +++ b/docs/figures/unweighted-sat/verification.svg @@ -0,0 +1 @@ + \ No newline at end of file diff --git a/docs/make.jl b/docs/make.jl index e14dd3a..363645a 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -42,6 +42,7 @@ makedocs(; "Rydberg Gadgets on Triangular Lattice" => "generated/trangular_Rydberg_example.md", "QUBO Gadgets on Triangular Lattice" => "generated/triangular_QUBO_example.md", ], + "Unweighted Search" => "unweighted_search.md", "Reference" => "ref.md", ], ) diff --git a/docs/src/unweighted_search.md b/docs/src/unweighted_search.md new file mode 100644 index 0000000..84b5002 --- /dev/null +++ b/docs/src/unweighted_search.md @@ -0,0 +1,118 @@ +# Four-pin unweighted lattice search + +`search_unweighted_gadgets` searches a finite lattice window for a four-pin +gadget whose reduced alpha tensor differs from the target by one constant. The +four pins and their outward directions must form the crossing geometry checked +by `check_crossing_frame`. + +The target graph is used only to obtain the target reduced alpha tensor and for +the unchanged final verifier. The search chooses the concrete pin locations, +pin directions, and lattice sites. `window_side` controls the finite search +window. The triangular lattice is the default; pass `Square()` to select KSG. + +The direct search enumerates the finite set of ordered pin locations and +outward lattice directions that satisfy the crossing geometry. For each layout, +atom count, and constant offset, a SAT instance chooses all occupied sites at +once. Its constraints enforce the target reduced alpha tensor and connectivity. +This does not require a known logical graph, a 28-site seed, or a historical +rewrite path. + +`min_vertices:max_vertices` is examined in increasing order, so reaching a +larger atom count means all scheduled smaller cases have already been rejected. +`max_evaluations` bounds SAT solver calls, not graph mutations. For each fixed +layout, atom count, and offset, satisfying site selections are blocked and the +solver is called again until that case is unsatisfiable. Every returned +candidate is independently checked by both `is_gadget_replacement` and +`check_crossing_frame`. + +`max_frame_evaluations` separately bounds pin-and-ray geometry checks. Increase +it together with `window_side` for exhaustive searches in larger windows. + +### Checkpoint and resume + +Pass `checkpoint_path` to persist deterministic frame, label-order, and offset +cursors after every `checkpoint_interval` SAT evaluations. Repeating the same +search call with the same path resumes from those cursors and restores its +evaluation and staged-filter counters. The target, lattice, vertex range, and +window size must match the saved search. + +`read_unweighted_search_checkpoint(path)` returns the saved `lattice`, +`atom_count`, `window_shape`, `frame_cursor`, `order_cursor`, `offset_cursor`, +`evaluated`, `frame_evaluated`, `frame_candidates`, `first_lower_rejected`, +`second_lower_rejected`, and `full_solves` fields for progress reporting. + +```julia +using GadgetSearch, Graphs + +target = SimpleGraph(4) +add_edge!(target, 1, 3) +add_edge!(target, 2, 4) + +result = search_unweighted_gadgets( + target, [1, 2, 3, 4], Triangular(); + min_vertices=4, + max_vertices=23, + max_evaluations=100_000, + max_frame_evaluations=10_000_000, + max_results=1, + window_side=8, +) +``` + +The keyword defaults are `min_vertices=5`, `max_vertices=13`, +`max_evaluations=2_000`, `max_frame_evaluations=1_000_000`, `max_results=1`, +and `window_side=4` for a four-pin target. The result contains the target, +ordered boundary, lattice name, accepted gadgets, number of SAT solver calls, and +one of four termination reasons: `:solution`, `:budget`, `:frame_budget`, or +`:search_space_exhausted`. + +## Joint occupancy-and-frame search + +`search_unweighted_gadget_joint` is the direct formulation used for blank-window +discovery. One CNF chooses the occupied sites, pins, and rays together. The +caller supplies one exact window shape, atom count, and tensor offset; the +function returns either a verifier-accepted `UnweightedGadget` or `nothing`. +`canonical_shift` and `first_ray` identify the symmetry-fixed first port for +this exact CNF instance. A large search should schedule those independent +instances externally with separate solver budgets. + +```julia +gadget = search_unweighted_gadget_joint( + target, [1, 2, 3, 4], Triangular(); + window_shape=(8, 8), atom_count=28, offset=10, + seconds=600, seed=1, canonical_shift=(0, 0), first_ray=1, +) +``` + +## Rewrite optimization + +`optimize_unweighted_gadget` is a local downstream atom-count optimizer. It uses +semantic rewrite rules rather than arbitrary vertex deletion: + +- contract a two-edge boundary tail while promoting its endpoint to the pin; +- contract opposite leaf pins (`P1`–`P3` or `P2`–`P4`) as one paired rewrite; +- move one frame pin by one lattice step and let fixed-frame SAT re-synthesize + every interior atom at a smaller atom count. + +Every accepted step passes the unchanged reduced-alpha verifier and all four +crossing-frame checks. Direct rules are explored to a closure rather than +greedily committing to the first smaller graph. Fixed-frame SAT is then tried +from the direct descendants, so a smaller direct dead end does not hide a +rewrite-and-resynthesize path. The result includes the best gadget reached in +the selected one-step frame neighborhood, a replayable before/after rewrite +trace, the number of fixed-frame SAT calls, and its termination reason. A fixed +point under these rules and budgets is not a global minimum certificate. Each +fixed-frame solve is capped by `max_sat_conflicts`; capped calls are counted in +`unresolved_sat_evaluations`. If the neighborhood is exhausted while any such +call remains unresolved, the termination reason is `:sat_unknown`, not +`:rewrite_fixed_point`. + +```julia +optimized = optimize_unweighted_gadget( + gadget, [1, 2, 3, 4]; + min_vertices=17, + max_sat_evaluations=256, + max_sat_conflicts=100_000, + host_radius=1, +) +``` diff --git a/docs/unweighted-sat-search-report.pdf b/docs/unweighted-sat-search-report.pdf new file mode 100644 index 0000000..149330c Binary files /dev/null and b/docs/unweighted-sat-search-report.pdf differ diff --git a/docs/unweighted-sat-search-report.typ b/docs/unweighted-sat-search-report.typ new file mode 100644 index 0000000..97d41d8 --- /dev/null +++ b/docs/unweighted-sat-search-report.typ @@ -0,0 +1,452 @@ +#set document( + title: "Certified Search and Rule-Based Optimization for Unweighted MIS Gadgets", + author: "GadgetSearch", +) +#set page( + paper: "a4", + margin: (left: 21mm, right: 21mm, top: 15mm, bottom: 17mm), + numbering: "1", + number-align: center + bottom, +) +#set text(font: "Libertinus Serif", size: 9.25pt) +#set math.equation(numbering: "(1)") +#set par(justify: true, leading: 0.58em) +#set heading(numbering: "1.1") +#set table(stroke: 0.45pt + rgb("#aeb7bf"), inset: (x: 5pt, y: 3.5pt)) +#show heading.where(level: 1): it => block( + above: 12pt, below: 5pt, + text(font: "Libertinus Sans", size: 14pt, weight: "semibold", it), +) +#show heading.where(level: 2): it => block( + above: 9pt, below: 3pt, + text(font: "Libertinus Sans", size: 11pt, weight: "semibold", it), +) +#show figure.caption: set text(size: 8.2pt) +#show raw: set text(font: "Libertinus Mono", size: 8pt) + +#let gray = rgb("#66727d") +#let fig(name, caption, width: 100%) = figure( + image("figures/unweighted-sat/" + name + ".svg", width: width), + caption: caption, +) + +#align(center)[ + #text(font: "Libertinus Serif Display", size: 19pt, weight: "semibold")[ + Certified Search and Rule-Based Optimization for Unweighted MIS Gadgets + ] + #v(2pt) + #text(size: 9pt)[GadgetSearch.jl] +] + +#v(4pt) +#block(inset: (x: 11mm, y: 7pt), stroke: ( + top: .6pt + gray, bottom: .6pt + gray, +))[ + #text(size: 9pt)[#text(font: "Libertinus Sans", weight: "semibold")[Abstract.] + We describe a certified synthesis framework for four-pin unweighted + maximum-independent-set gadgets on the triangular and king's graph + lattices. A joint SAT formulation discovers a first gadget by choosing its + occupied sites and four ports in one finite window. A complementary + frame-first formulation enumerates port geometry and solves a smaller + occupancy problem for each frame. Explicit independent-set witnesses impose + lower bounds on the reduced alpha tensor, frontier binary decision diagrams + impose upper bounds, and direct layered reachability imposes connectivity. + The second stage treats successful reductions as evidence from which + to extract reusable, interface-preserving rewrite schemas. Their common + invariant is equality, up to one constant, of the conditioned MIS signature + seen across a small separator. Exact local replacements and + interface-constrained SAT re-synthesis are instances of this rule. Every + application is certified by the unchanged verifier. A reduction from a + blank-window 28-site CROSS to 23 sites serves as a regression case, not as + the definition of the optimization algorithm.] +] + += Problem formulation + +Let $R$ be a target graph with four ordered boundary vertices +$partial R=(b_1,b_2,b_3,b_4)$. For a boundary state +$sigma=(sigma_1,...,sigma_4) in {0,1}^4$, the value $sigma_i=1$ requires +$b_i$ to belong to the independent set. Define +$ + alpha(R)_sigma = max{abs(I): I " independent in " R, + I ∩ partial R = {b_i:sigma_i=1}}. +$ +An infeasible boundary state has value $-infinity$. The *reduced alpha tensor* +$tilde(alpha)(R)$ removes entries dominated by a proper subconfiguration. +Denote the target tensor by $T=tilde(alpha)(R)$. + +Let $Lambda$ be a lattice and let $H$ be the blockade graph of a finite lattice +window. A selected site set $S subset.eq V(H)$ induces the replacement graph +$G=H[S]$. Four ordered selected sites $P=(p_1,p_2,p_3,p_4)$ represent the +boundary vertices, and lattice directions $bold(r)=(r_1,r_2,r_3,r_4)$ describe +the exterior wires. The search seeks $S$, $P$, $bold(r)$, and an integer offset +$c$ satisfying +$ + tilde(alpha)(G) = T+c. +$ + +For each pin, let $s_i=p_i+r_i$ be the first interface site outside the gadget. +The frame must satisfy four geometric conditions: + +- *(G1)* every $s_i$ is a strict vertex of the convex hull of + $S union {s_1,s_2,s_3,s_4}$; +- *(G2)* the cyclic hull order alternates the interface pairs $(s_1,s_3)$ and + $(s_2,s_4)$; +- *(G3)* every $r_i$ points strictly outward from the interface centroid; +- *(G4)* the exterior rays are unobstructed, touch only their own pins, and + remain pairwise outside blockade range. + +The selected graph must contain exactly $N$ sites, contain all four pins, and +be connected. A candidate is accepted only when the unchanged +`is_gadget_replacement` routine independently recomputes the tensor, agrees on +$c$, and the geometry checker confirms G1--G4. + += Two-stage algorithm + +The recommended workflow has two stages: *search* first produces one or more +certified seed gadgets, and *optimization* then attempts to reduce each seed. +The public API exposes the two stages separately so that searches and optimizer +budgets can be scheduled independently; a complete synthesis run calls them in +this order. + +#fig("pipeline", [The complete algorithm is Stage I search followed by Stage II +optimization. Joint and frame-first SAT are alternative engines inside the +search stage, not alternatives to optimization.]) + +Stage I has two search engines that share the tensor and connectivity encodings: + +- *Joint occupancy-and-frame SAT* chooses $S$, $P$, and $bold(r)$ in one CNF. + It is the primary blank-window discovery method. +- *Frame-first SAT* enumerates $(P,bold(r))$ geometrically and lets SAT choose + only $S$. It is suited to systematic finite-window traversal and local + re-synthesis around a known gadget. + +Stage II consumes every verifier-valid `UnweightedGadget` returned by Stage I. +It never assumes that the starting gadget came from a particular solver run. +If no smaller certified descendant exists within the rule and SAT budgets, the +seed itself is the final result of the two-stage algorithm. + += Joint occupancy-and-frame SAT + +Fix a lattice window $W$, atom count $N$, and offset $c$. For every site +$v in W$, introduce a selection variable $x_v$. For every geometrically +eligible site-direction pair and boundary label $i$, introduce a frame variable +$f_(i,v,d)$. The joint CNF imposes: + +1. $sum_v x_v=N$; +2. exactly one frame choice for each label $i$; +3. each chosen pin is selected and the four pins are distinct; +4. local corridor, outward-direction, and alternating-interface constraints; +5. the reduced-alpha tensor constraints described below; +6. connectivity of the selected induced graph. + +Only eligible port choices receive variables. This sparse allocation is much +smaller than allocating every $(i,v,d)$ combination and forbidding most of them. + +One reference port is fixed by a canonical translation and ray direction. This +breaks translation and rotation symmetry, but it also defines the represented +finite instance. A known abstract gadget is absent if no symmetry-equivalent +embedding fits the anchored window. Canonical placement is therefore part of +the instance specification, not an innocuous implementation detail. + +The formula is emitted as DIMACS and solved by Kissat. A model is materialized +and checked geometrically before the reduced tensor is recomputed. A rejected +model contributes a blocker and the solver is restarted; no rejected candidate +is reported as a gadget. + +== Solver portfolios + +Blank-window joint instances can be highly seed-sensitive even when they +represent the same finite search space. Independent Kissat seeds are therefore +a useful source of parallelism. A timed-out run is recorded as `UNKNOWN` and +may be rescheduled; it is never interpreted as an UNSAT certificate. + += Frame-first search + +Frame-first search moves geometry outside the occupancy CNF. It constructs a +compatibility graph whose vertices are locally valid $(p,r)$ candidates. Four +mutually compatible candidates form a possible four-port frame. Valid label +orders are then checked against the complete G1--G4 conditions and canonicalized +under the lattice symmetries. + +#fig("frame-enumeration", [Frame-first construction. Compatible port candidates +form a four-clique; a valid ordered frame determines the allowed host before +occupancy variables are introduced.]) + +For a surviving frame, each lattice site is tested against its pin corridors +and hull constraints. The allowed set $U(P,bold(r))$ defines the host +$H=G[U]$. The fixed-frame CNF contains site-selection, tensor, and connectivity +variables, but no frame-choice variables. + +For the known 23-site triangular CROSS frame, Kissat selects the occupied sites +and returns a verifier-valid gadget at offset 7. The hard part of a blank +frame-first search is reaching this rare frame, not solving its occupancy +instance. + +The production frame traversal has deterministic cursors and periodically +writes checkpoints containing the current frame, label order, offset, and +filter counters. Repeating the same bounded search with that checkpoint resumes +the traversal. This makes long finite traversals resumable, although their +total geometry space can still be very large. + += SAT encoding for a fixed host + +Fix $(H,P,N,c)$. The formula is satisfiable exactly when $H$ contains a +connected $N$-site induced subgraph with reduced alpha tensor $T+c$. + +== Atom selection + +For every host vertex $v$, introduce $x_v$, with +$ + x_(p_i)=1 quad (i=1,...,4), + quad sum_(v in V(H)) x_v=N. +$ + +== Direct layered connectivity + +Let $z_v^t$ mean that selected vertex $v$ is reachable from $p_1$ by a path of +at most $t$ host edges. The initial layer contains only the root, +$ + z_v^0 arrow.l.r.double v=p_1. +$ +For $t=0,...,N-2$, direct propagation is +$ + z_v^(t+1) arrow.l.r.double + (z_v^t or (x_v and or.big_(u in N_H(v)) z_u^t)). +$ +Finally, +$ + x_v => z_v^(N-1) quad (v in V(H)). +$ +No per-edge arrival variables are used. For a host with $h=abs(V(H))$ vertices +and $m=abs(E(H))$ undirected edges, this encoding contributes $N h$ auxiliary +variables and +$ + 2h+(N-1)(3h+2m) +$ +clauses. The dependence on host edges is linear, matching the sparse TLSG and +KSG hosts. + +#fig("connectivity-layers", [Direct layered reachability for a selected path +$p_1-u-v$ and an isolated host vertex $w$. Reachability persists and propagates +through selected vertices; $w$ cannot satisfy the final implication.], width: 78%) + +The depth $N-1$ is exact: every vertex in a connected $N$-vertex graph has a +simple path of length at most $N-1$ from $p_1$, whereas a selected vertex in +another component can never become reachable. + +== From the reduced tensor to bounds + +For a selected set $S$ and boundary state $sigma$, let $A_S(sigma)$ be the +largest independent-set size with exactly the pins indicated by $sigma$. +Define the monotone completion +$ + C_T(sigma)=max{T_tau: tau subset.eq sigma, T_tau != -infinity}. +$ +The equality $tilde(alpha)(G[S])=T+c$ is imposed by +$ + A_S(sigma) >= T_sigma+c quad "for finite " T_sigma, +$ +and +$ + A_S(sigma) <= C_T(sigma)+c quad "for every " sigma. +$ +Only lower and upper states not implied by other boundary states are emitted. + +== Lower bounds: explicit witnesses + +For each essential finite $T_sigma$, introduce witness variables +$y_v^sigma$. The clauses impose +$ + y_v^sigma => x_v, + quad y_u^sigma+y_v^sigma <= 1 quad ({u,v} in E(H)), +$ +$ + y_(p_i)^sigma=sigma_i, + quad sum_(v in V(H)) y_v^sigma=T_sigma+c. +$ +Thus the selected graph explicitly contains an independent set of the required +size and boundary state. + +== Upper bounds: frontier BDDs + +Upper bounds must exclude every oversized independent set. The implementation +constructs a frontier binary decision diagram under a fixed vertex order. A +path either skips or takes the current host vertex. A take edge is enabled only +when the site is selected and no taken frontier vertex conflicts with it. + +After each prefix, BDD states with the same active frontier are merged. The +count is capped at $C_T(sigma)+c+1$, and every terminal state at the cap is +forbidden. The implementation tries three lattice projections and their +reversals and retains the smallest of the six BDDs. + +#fig("bdd-example", [Frontier BDD for the path $a-b-c$ with upper bound one. +Skip and take branches represent every independent set; the two-vertex terminal +is forbidden.], width: 72%) + +For fixed $(H,P,N,c)$, the complete occupancy formula is +$ + Phi = Phi_"select" and Phi_"connect" + and ∧_(sigma in L) Phi_"witness"^sigma + and ∧_(sigma in U) Phi_"BDD"^sigma, +$ +where $L$ and $U$ are the essential lower and upper boundary states. + += Independent certification + +A satisfying assignment is not returned directly. The selected coordinates +are materialized as the exact blockade graph, after which the production +verifier recomputes $tilde(alpha)(G)$ and checks the common offset. Connectivity +and G1--G4 are checked independently. + +#fig("verification", [The SAT solver proposes; the unchanged verifier accepts. +The same certification boundary is used during discovery and after every +rewrite.]) + +This separation protects the result from an incomplete SAT-side optimization: +auxiliary encodings may change, but the semantic acceptance test does not. + += Reusable rewrite principle and implemented optimizer + +A successful reduction should not be stored merely as a sequence of edits. +The useful object is the *semantic reason* why the changed region can be +replaced in any compatible context. For MIS gadgets, that reason is a boundary +signature across a small separator. + +== Interface signatures + +Let $Q$ be a connected subgraph of a certified gadget and let +$C=(c_1,...,c_k)$ contain every vertex through which $Q$ meets the unchanged +exterior. For an interface state $tau in {0,1}^k$, define +$ + A_Q^C(tau)=max{abs(I): I " independent in " Q, + I ∩ C={c_i:tau_i=1}}. +$ +The vector $"Sig"_C(Q)=(A_Q^C(tau))_tau$ is everything the exterior needs to +know about $Q$ when maximizing an independent set. A local replacement +$Q arrow.r Q'$ is context-independent when +$ + A_(Q')^C(tau)=A_Q^C(tau)+delta + quad "for every feasible " tau, +$ +and the infeasible-state pattern is unchanged. Gluing either patch to the same +exterior then changes every global conditioned optimum by the same $delta$. +Consequently, the reduced alpha tensor of the whole gadget changes only by the +allowed constant offset. + +The signature condition is necessary but not sufficient for a lattice rewrite. +The replacement must also preserve the attachment coordinates or port roles, +introduce no unintended blockade edge to the exterior, maintain connectivity, +and satisfy G1--G4. The unchanged whole-gadget verifier remains the final +certificate. + +== Generalizing one successful case + +A successful before/after pair suggests a reusable rule only when its changed +region can be isolated behind a small separator and the two local interface +signatures differ by one constant. Absolute coordinates must then be replaced +by structural premises such as interface order, adjacency, pin role, channel +pairing, and lattice direction. This is the criterion for adding future rules: +the proof obligation is the interface contract, not reproduction of one edit +sequence. The current code does not automatically mine or minimize schemas. + +== Operators implemented today + +The current optimizer implements two operator families consistent with the +signature principle. + +- *Exact local replacement.* Two closed-form rules contract an even boundary + tail or a valid opposite pair of leaf pins. Each proposed whole gadget is + certified; the implementation does not yet enumerate arbitrary equivalent + local patches from their signatures. +- *Interface-constrained re-synthesis.* Hold the separator, external port roles, + and required signature fixed, then ask SAT to synthesize a smaller interior. + Moving a pin by one lattice step is one proposal mechanism for changing the + interface geometry, not the semantic rewrite rule. The rule is to preserve + the interface contract while re-solving the interior. + +The optimizer searches the graph generated by these operators. It explores +all certified direct descendants to a closure and only then spends SAT budget +on re-synthesis. Exact coordinate, pin, and ray state keys merge identical +converging paths. This non-greedy search is essential: the smallest immediate +child may be a dead end, while a different equivalent child exposes a stronger +subsequent rewrite. + +== What the CROSS case teaches + +The reduction from the independently discovered 28-site CROSS to a 23-site +gadget exercises both implemented operator families: certified structural +contractions and one interface-constrained interior re-synthesis. Its general +lesson is to preserve an interface contract, explore alternative certified +descendants, and invoke exact re-synthesis when a structural step changes the +available interior. The particular atom counts and rule names are regression +evidence, not assumptions of a new search. + +The same second stage is applied to CROSS+EDGE. Stage I finds a verified +9-site gadget; Stage II exhausts the direct contraction closure and a budget of +256 fixed-frame SAT evaluations without finding a smaller certified descendant. +The reported 9-site result is therefore the best result under the current rules +and budget, not a search-only result and not a proof of global minimality. + +The final result is certified relative to the available rule schemas, host +radius, minimum atom count, and SAT budget. It is not a proof that no smaller +gadget exists. + += Results and regression coverage + +#figure( + table( + columns: (1.35fr, 1.35fr, 1.35fr, 1.35fr, 2.3fr), + align: (left, center, center, center, left), + fill: (x, y) => if y == 0 { rgb("#edf4fa") } else { white }, + table.header( + [*target*], [*Stage I search*], [*Stage II optimize*], + [*offset*], [*final status*], + ), + [CROSS], [28 sites], [23 sites], [$10 arrow.r 7$], + [reduced; verified], + [CROSS+EDGE], [9 sites], [9 sites], [$1 arrow.r 1$], + [unchanged after budget; verified], + ), + caption: [End-to-end results of the complete two-stage algorithm.], +) + +Correctness is exercised at several levels: + +- direct layered connectivity is compared with explicit connected-subset + enumeration on small host graphs; +- fixed-frame and joint positive controls pass the unchanged verifier; +- Stage I search covers both CROSS and CROSS+EDGE, and Stage II optimization is + run on both results; +- the historical 28-site CROSS satisfies the production CNF and verifier; +- a regression starts from that real 28-site result and checks every certified + rewrite application down to 23 sites; +- the focused unweighted-search suite and the complete package test suite pass. + += Completeness and practical limits + +The procedure is complete only within its explicitly scheduled finite space. + +- A joint instance fixes one window, $N$, $c$, canonical anchor, and first ray. +- A frame-first run fixes a finite window, atom-count range, frame budget, and + SAT-call budget. +- Kissat `UNKNOWN` records expired work, not a negative proof. +- The rewrite optimizer explores a bounded signature-preserving rule + neighborhood and does not + certify global minimality. + +Frame-first enumeration can encounter millions of easy UNSAT frames before a +useful frame, with a small hard tail dominating wall time. Joint SAT avoids that +explicit traversal but exposes one larger, seed-sensitive CNF. The two modes are +therefore complementary rather than competing implementations. + += Recommended workflow + +For a new four-pin unweighted target: + +1. *Stage I -- search:* compute the target tensor, choose TLSG or KSG, and run + joint or frame-first SAT until at least one seed gadget is independently + verified; +2. *Stage II -- optimize:* run exact local replacements and + interface-constrained re-synthesis on every distinct seed, retaining the + smallest verified descendant or the unchanged seed at a fixed point; diff --git a/src/GadgetSearch.jl b/src/GadgetSearch.jl index 430f671..1269d70 100644 --- a/src/GadgetSearch.jl +++ b/src/GadgetSearch.jl @@ -9,6 +9,7 @@ using IterTools using Combinatorics using Random using Serialization +import Kissat_jll # Tensor network dependencies for alpha tensor computation using GenericTensorNetworks: GenericTensorNetwork, IndependentSet, SizeMax, solve @@ -26,6 +27,7 @@ include("graphio/savegraph.jl") include("graphio/udg.jl") include("utils/ruleio.jl") include("core/unweighted_search.jl") +include("core/unweighted_sat.jl") include("core/search.jl") include("utils/gadget.jl") include("utils/visualize.jl") @@ -72,6 +74,13 @@ export is_gadget_replacement # Unweighted search export UnweightedGadget +export UnweightedSearchResult +export UnweightedRewriteStep +export UnweightedOptimizationResult export search_unweighted_gadgets +export search_unweighted_gadget_joint +export optimize_unweighted_gadget +export read_unweighted_search_checkpoint +export check_crossing_frame end # module diff --git a/src/core/unweighted_sat.jl b/src/core/unweighted_sat.jl new file mode 100644 index 0000000..3f77c2f --- /dev/null +++ b/src/core/unweighted_sat.jl @@ -0,0 +1,1238 @@ +mutable struct _SatCnf + variables::Int + clauses::Vector{Vector{Int}} +end + +struct _SatFrameContext + frame + coordinates::Vector{_LatticeCoordinate} + boundary::Vector{Int} + edge_list::Vector{Tuple{Int, Int}} + order::Vector{Int} + layers::Vector{Vector{Tuple{Int, Int, Bool}}} + adjacent::Vector{Vector{Int}} +end + +mutable struct _SatSearchStats + frame_candidates::Int + first_lower_rejected::Int + second_lower_rejected::Int + full_solves::Int +end + +struct _SatSearchCheckpoint + target::Vector{Float64} + lattice::Symbol + min_vertices::Int + max_vertices::Int + window_side::Int + atom_count::Int + window_shape::Tuple{Int, Int} + frame_cursor::Int + order_cursor::Int + offset_cursor::Int + evaluated::Int + frame_evaluated::Int + stats::_SatSearchStats +end + +_SatSearchStats() = _SatSearchStats(0, 0, 0, 0) + +function _write_sat_search_checkpoint(path, checkpoint) + open(path, "w") do stream + serialize(stream, checkpoint) + end +end + +function _read_sat_search_checkpoint(path) + open(path) do stream + return deserialize(stream) + end +end + +"""Read progress and staged-filter counters from a saved unweighted search.""" +function read_unweighted_search_checkpoint(path::String) + checkpoint = _read_sat_search_checkpoint(path) + return ( + lattice=checkpoint.lattice, + atom_count=checkpoint.atom_count, + window_shape=checkpoint.window_shape, + frame_cursor=checkpoint.frame_cursor, + order_cursor=checkpoint.order_cursor, + offset_cursor=checkpoint.offset_cursor, + evaluated=checkpoint.evaluated, + frame_evaluated=checkpoint.frame_evaluated, + frame_candidates=checkpoint.stats.frame_candidates, + first_lower_rejected=checkpoint.stats.first_lower_rejected, + second_lower_rejected=checkpoint.stats.second_lower_rejected, + full_solves=checkpoint.stats.full_solves, + ) +end + +_SatCnf() = _SatCnf(0, Vector{Int}[]) + +function _sat_variable!(cnf::_SatCnf) + cnf.variables += 1 + return cnf.variables +end + +_sat_clause!(cnf::_SatCnf, literals::Int...) = push!(cnf.clauses, collect(literals)) + +function _sat_at_most!(cnf::_SatCnf, literals::Vector{Int}, bound::Int) + bound >= length(literals) && return + bound >= 0 || return _sat_clause!(cnf) + if bound == 0 + foreach(literal -> _sat_clause!(cnf, -literal), literals) + return + end + counters = [ + [_sat_variable!(cnf) for _ in 1:bound] + for _ in 1:length(literals)-1 + ] + for row in 1:length(literals)-1 + _sat_clause!(cnf, -literals[row], counters[row][1]) + end + for row in 2:length(literals)-1 + _sat_clause!(cnf, -counters[row-1][1], counters[row][1]) + end + for column in 2:bound + for row in column:length(literals)-1 + _sat_clause!( + cnf, -literals[row], -counters[row-1][column-1], + counters[row][column], + ) + end + for row in column+1:length(literals)-1 + _sat_clause!(cnf, -counters[row-1][column], counters[row][column]) + end + end + for row in bound+1:length(literals) + _sat_clause!(cnf, -literals[row], -counters[row-1][bound]) + end +end + +function _sat_exactly!(cnf::_SatCnf, literals::Vector{Int}, count::Int) + _sat_at_most!(cnf, literals, count) + _sat_at_most!(cnf, -literals, length(literals) - count) +end + +function _sat_exactly_one!(cnf::_SatCnf, literals::Vector{Int}) + _sat_clause!(cnf, literals...) + _sat_at_most!(cnf, literals, 1) +end + +function _independent_set_bdd(coordinates, edge_list) + adjacent = [Set{Int}() for _ in coordinates] + for (first, second) in edge_list + push!(adjacent[first], second) + push!(adjacent[second], first) + end + canonical = _offset_to_axial.(coordinates) + orderings = Vector{Vector{Int}}() + for projection in ( + point -> (point[1], point[2]), + point -> (point[2], point[1]), + point -> (point[1] + point[2], point[1]), + ) + order = sortperm(canonical; by=projection) + push!(orderings, order, reverse(order)) + end + function build_layers(order) + position = zeros(Int, length(coordinates)) + for (step, vertex) in enumerate(order) + position[vertex] = step + end + last_neighbor = [maximum([position[vertex]; position[collect(adjacent[vertex])]]) + for vertex in eachindex(coordinates)] + layer = [Int[]] + layers = Vector{Vector{Tuple{Int, Int, Bool}}}() + layer_sizes = [1] + for (step, vertex) in enumerate(order) + next_layer = Vector{Vector{Int}}() + next_index = Dict{Tuple{Vararg{Int}}, Int}() + arcs = Tuple{Int, Int, Bool}[] + for (source, occupied) in enumerate(layer), take in (false, true) + take && any(neighbor -> neighbor in occupied, adjacent[vertex]) && continue + next_occupied = take ? [occupied; vertex] : occupied + frontier = Tuple(sort!(filter(v -> last_neighbor[v] > step, next_occupied))) + destination = get!(next_index, frontier) do + push!(next_layer, collect(frontier)) + length(next_layer) + end + push!(arcs, (source, destination, take)) + end + push!(layers, arcs) + push!(layer_sizes, length(next_layer)) + layer = next_layer + end + return sum(layer_sizes), order, layers + end + best = build_layers(first(orderings)) + for order in Iterators.drop(orderings, 1) + candidate = build_layers(order) + candidate[1] < best[1] && (best = candidate) + end + _, order, layers = best + return order, layers +end + +function _target_completion(target) + return [maximum(target[subset+1] for subset in 0:state + if subset & ~state == 0 && isfinite(target[subset+1])) + for state in 0:length(target)-1] +end + +function _essential_lower_states(target) + states = Int[] + for state in 0:length(target)-1 + value = target[state+1] + isfinite(value) || continue + implied = [target[superset+1] - count_ones(superset ⊻ state) + for superset in state+1:length(target)-1 + if state & ~superset == 0 && isfinite(target[superset+1])] + (isempty(implied) || value > maximum(implied)) && push!(states, state) + end + return states +end + +function _essential_upper_states(completion) + states = Int[] + for state in 0:length(completion)-1 + implied = [completion[subset+1] + count_ones(state ⊻ subset) + for subset in 0:state-1 if subset & ~state == 0] + (isempty(implied) || completion[state+1] < minimum(implied)) && + push!(states, state) + end + return states +end + +function _centered_offsets(offsets, center=(first(offsets) + last(offsets)) / 2) + return sort!(collect(offsets); by=offset -> (abs(offset - center), offset)) +end + +function _add_independent_set_upper_bound!( + cnf, selected, pins, order, layers, boundary_state, bound, +) + root = _sat_variable!(cnf) + _sat_clause!(cnf, root) + reachable = Dict((1, 0) => root) + pin_slots = Dict(pin => slot - 1 for (slot, pin) in enumerate(pins)) + for (step, (vertex, arcs)) in enumerate(zip(order, layers)) + next_reachable = Dict{Tuple{Int, Int}, Int}() + required = haskey(pin_slots, vertex) ? + !iszero(boundary_state & (1 << pin_slots[vertex])) : nothing + for (source, destination, take) in arcs + required !== nothing && take != required && continue + for ((source_state, count), source_variable) in reachable + source_state == source || continue + next_count = min(bound + 1, count + Int(take)) + destination_variable = get!(next_reachable, (destination, next_count)) do + _sat_variable!(cnf) + end + take ? _sat_clause!( + cnf, -source_variable, -selected[vertex], destination_variable, + ) : _sat_clause!(cnf, -source_variable, destination_variable) + end + end + reachable = next_reachable + end + for ((_, count), variable) in reachable + count == bound + 1 && _sat_clause!(cnf, -variable) + end +end + +function _add_selected_connectivity!( + cnf, selected, adjacent, root, selected_count, +) + reachable = [_sat_variable!(cnf) for _ in selected] + for vertex in eachindex(reachable) + _sat_clause!(cnf, vertex == root ? reachable[vertex] : -reachable[vertex]) + end + for _ in 1:selected_count-1 + next_reachable = [_sat_variable!(cnf) for _ in selected] + for vertex in eachindex(selected) + _sat_clause!(cnf, -reachable[vertex], next_reachable[vertex]) + _sat_clause!(cnf, -next_reachable[vertex], selected[vertex]) + _sat_clause!( + cnf, -next_reachable[vertex], reachable[vertex], + (reachable[neighbor] for neighbor in adjacent[vertex])..., + ) + for neighbor in adjacent[vertex] + _sat_clause!( + cnf, -selected[vertex], -reachable[neighbor], + next_reachable[vertex], + ) + end + end + reachable = next_reachable + end + for vertex in eachindex(selected) + _sat_clause!(cnf, -selected[vertex], reachable[vertex]) + end +end + +function _prepare_sat_frame(lattice, frame) + coordinates = sort(copy(frame.allowed)) + patch = _LatticePatch(coordinates, copy(frame.pins), copy(frame.rays)) + host, boundary, _ = _materialize_lattice_patch(lattice, patch) + edge_list = [(src(edge), dst(edge)) for edge in edges(host)] + order, layers = _independent_set_bdd(coordinates, edge_list) + adjacent = [collect(neighbors(host, vertex)) for vertex in vertices(host)] + return _SatFrameContext( + frame, coordinates, boundary, edge_list, order, layers, adjacent, + ) +end + +function _solve_fixed_crossing_sat(target, lattice, frame, atom_count, offset) + return _solve_fixed_crossing_sat( + target, lattice, _prepare_sat_frame(lattice, frame), atom_count, offset, + ) +end + +function _solve_fixed_crossing_sat( + target, lattice, context::_SatFrameContext, atom_count, offset, +) + solver, selected = _fixed_crossing_sat_problem( + target, context, atom_count, offset, + ) + return _solve_next_fixed_crossing_sat!( + solver, selected, target, lattice, context, + ) +end + +function _fixed_crossing_sat_problem(target, context, atom_count, offset) + cnf, selected = _fixed_crossing_sat_cnf( + target, context, atom_count, offset, + ) + return _new_sat_solver(cnf), selected +end + +function _fixed_crossing_sat_cnf(target, context, atom_count, offset) + coordinates = context.coordinates + pins = context.boundary + cnf = _SatCnf() + selected = [_sat_variable!(cnf) for _ in coordinates] + foreach(pin -> _sat_clause!(cnf, selected[pin]), pins) + _sat_exactly!(cnf, selected, atom_count) + + for state in _essential_lower_states(target) + _add_lower_state_constraint!( + cnf, selected, target, context, state, offset, + ) + end + + completion = _target_completion(target) + for state in _essential_upper_states(completion) + _add_independent_set_upper_bound!( + cnf, selected, pins, context.order, context.layers, state, + Int(completion[state+1] + offset), + ) + end + _add_selected_connectivity!( + cnf, selected, context.adjacent, pins[1], atom_count, + ) + return cnf, selected +end + +function _add_lower_state_constraint!( + cnf, selected, target, context, state, offset, +) + witness = [_sat_variable!(cnf) for _ in context.coordinates] + for vertex in eachindex(context.coordinates) + _sat_clause!(cnf, -witness[vertex], selected[vertex]) + end + for (slot, pin) in enumerate(context.boundary) + _sat_clause!( + cnf, iszero(state & (1 << (slot - 1))) ? + -witness[pin] : witness[pin], + ) + end + for (first, second) in context.edge_list + _sat_clause!(cnf, -witness[first], -witness[second]) + end + _sat_exactly!(cnf, witness, Int(target[state+1] + offset)) +end + +function _add_joint_lower_state_constraint!( + cnf, selected, target, edge_list, frame_choices, state, offset, +) + witness = [_sat_variable!(cnf) for _ in selected] + for vertex in eachindex(selected) + _sat_clause!(cnf, -witness[vertex], selected[vertex]) + end + for label in eachindex(frame_choices), + (vertex, _, port) in frame_choices[label] + occupied = !iszero(state & (1 << (label - 1))) + _sat_clause!(cnf, -port, + occupied ? witness[vertex] : -witness[vertex]) + end + for (first, second) in edge_list + _sat_clause!(cnf, -witness[first], -witness[second]) + end + _sat_exactly!(cnf, witness, Int(target[state+1] + offset)) +end + +function _add_joint_independent_set_upper_bound!( + cnf, selected, choices_by_slot_vertex, order, layers, boundary_state, bound, +) + maximum_remaining = [Dict{Int, Int}() for _ in 1:length(order)+1] + maximum_remaining[end][1] = 0 + for step in length(order):-1:1 + for (source, destination, take) in layers[step] + haskey(maximum_remaining[step+1], destination) || continue + value = Int(take) + maximum_remaining[step+1][destination] + maximum_remaining[step][source] = max( + get(maximum_remaining[step], source, -1), value, + ) + end + end + root = _sat_variable!(cnf) + _sat_clause!(cnf, root) + reachable = [(1, 0, root)] + for (step, (vertex, arcs)) in enumerate(zip(order, layers)) + required_taken = Int[] + required_skipped = Int[] + for label in eachindex(choices_by_slot_vertex) + destination = iszero(boundary_state & (1 << (label - 1))) ? + required_skipped : required_taken + append!(destination, choices_by_slot_vertex[label][vertex]) + end + next_reachable = Tuple{Int, Int, Int}[] + next_index = Dict{Tuple{Int, Int}, Int}() + for (source, destination, take) in arcs + for (source_state, count, source_variable) in reachable + source_state == source || continue + next_count = min(bound + 1, count + Int(take)) + next_count + get(maximum_remaining[step+1], destination, -1) < + bound + 1 && continue + key = (destination, next_count) + destination_variable = get(next_index, key, 0) + if iszero(destination_variable) + destination_variable = _sat_variable!(cnf) + next_index[key] = destination_variable + push!(next_reachable, (destination, next_count, destination_variable)) + end + if take + _sat_clause!( + cnf, -source_variable, -selected[vertex], + required_skipped..., destination_variable, + ) + else + _sat_clause!( + cnf, -source_variable, required_taken..., + destination_variable, + ) + end + end + end + reachable = next_reachable + end + for (_, count, variable) in reachable + count == bound + 1 && _sat_clause!(cnf, -variable) + end +end + +function _joint_frame_choice_is_outward(lattice, interface, direction, interfaces) + pin = (interface[1] - direction[1], interface[2] - direction[2]) + geometry = _geometry_coordinate.(Ref(lattice), _from_canonical.(Ref(lattice), interfaces)) + pin_geometry = _geometry_coordinate(lattice, _from_canonical(lattice, pin)) + hull = _strict_convex_hull([geometry; pin_geometry]) + all(in(hull), geometry) || return false + sum_q = sum(first, interfaces) + sum_r = sum(last, interfaces) + out_q = 4interface[1] - sum_q + out_r = 4interface[2] - sum_r + if lattice isa Square + return out_q * direction[1] + out_r * direction[2] > 0 + end + out_x = 2out_q + out_r + direction_x = 2direction[1] + direction[2] + return out_x * direction_x + 3out_r * direction[2] > 0 +end + +function _sat_group_indicator!(cnf, literals) + indicator = _sat_variable!(cnf) + for literal in literals + _sat_clause!(cnf, -literal, indicator) + end + _sat_clause!(cnf, -indicator, literals...) + return indicator +end + +function _interfaces_form_alternating_quadrilateral(interfaces) + length(unique(interfaces)) == 4 || return false + return ( + _orientation(interfaces[1], interfaces[3], interfaces[2]) * + _orientation(interfaces[1], interfaces[3], interfaces[4]) < 0 && + _orientation(interfaces[2], interfaces[4], interfaces[1]) * + _orientation(interfaces[2], interfaces[4], interfaces[3]) < 0 + ) +end + +function _joint_rays_touch(start1, direction1, start2, direction2, directions) + return any([_LatticeCoordinate[(0, 0)]; directions]) do offset + _rays_touch(start1, direction1, start2, direction2, offset) + end +end + +function _joint_site_blocks_ray(site, pin, direction, directions) + site == pin && return false + start = (pin[1] + direction[1], pin[2] + direction[2]) + return any([_LatticeCoordinate[(0, 0)]; directions]) do offset + _point_on_ray( + (site[1] + offset[1], site[2] + offset[2]), start, direction, + ) + end +end + +function _add_joint_frame_choices!( + cnf, selected, coordinates, lattice, first_direction_index, +) + directions = _lattice_directions(lattice) + canonical = _canonical_coordinate.(Ref(lattice), coordinates) + origin = only(findall(==((0, 0)), canonical)) + frame_choices = [Tuple{Int, Int, Int}[] for _ in 1:4] + push!(frame_choices[1], ( + origin, first_direction_index, _sat_variable!(cnf), + )) + fixed_start = directions[first_direction_index] + for label in 2:4, vertex in eachindex(coordinates), + direction_index in eachindex(directions) + direction = directions[direction_index] + pin = canonical[vertex] + start = (pin[1] + direction[1], pin[2] + direction[2]) + (vertex == origin || _joint_rays_touch( + fixed_start, directions[first_direction_index], start, direction, + directions, + )) && continue + push!(frame_choices[label], + (vertex, direction_index, _sat_variable!(cnf))) + end + + choices_by_slot_vertex = [ + [Int[] for _ in coordinates] for _ in 1:4 + ] + for label in eachindex(frame_choices) + _sat_exactly_one!(cnf, last.(frame_choices[label])) + for (vertex, _, variable) in frame_choices[label] + push!(choices_by_slot_vertex[label][vertex], variable) + _sat_clause!(cnf, -variable, selected[vertex]) + end + end + for first_label in 1:3, second_label in first_label+1:4, + vertex in eachindex(coordinates), + first in choices_by_slot_vertex[first_label][vertex], + second in choices_by_slot_vertex[second_label][vertex] + _sat_clause!(cnf, -first, -second) + end + + interface_variables = [Dict{_LatticeCoordinate, Int}() for _ in 1:4] + interface_order = [_LatticeCoordinate[] for _ in 1:4] + records = [ + Dict{_LatticeCoordinate, Vector{Tuple{Int, Int}}}() for _ in 1:4 + ] + for label in eachindex(frame_choices) + for (vertex, direction_index, choice) in frame_choices[label] + pin = canonical[vertex] + direction = directions[direction_index] + interface = (pin[1] + direction[1], pin[2] + direction[2]) + if !haskey(records[label], interface) + records[label][interface] = Tuple{Int, Int}[] + push!(interface_order[label], interface) + end + push!(records[label][interface], (direction_index, choice)) + end + for interface in interface_order[label] + choices = last.(records[label][interface]) + variable = _sat_group_indicator!(cnf, choices) + interface_variables[label][interface] = variable + end + end + + first_interface = only(interface_order[1]) + adjacent_pairs = Tuple{_LatticeCoordinate, _LatticeCoordinate, Int}[] + for adjacent_first in interface_order[2], adjacent_second in interface_order[4] + adjacent_first < adjacent_second || continue + pair = _sat_variable!(cnf) + first_variable = interface_variables[2][adjacent_first] + second_variable = interface_variables[4][adjacent_second] + _sat_clause!(cnf, -pair, first_variable) + _sat_clause!(cnf, -pair, second_variable) + _sat_clause!(cnf, -first_variable, -second_variable, pair) + push!(adjacent_pairs, (adjacent_first, adjacent_second, pair)) + end + for opposite in interface_order[3] + opposite_variable = interface_variables[3][opposite] + allowed_pairs = Int[] + for (adjacent_first, adjacent_second, pair) in adjacent_pairs + interfaces = [ + first_interface, adjacent_first, opposite, adjacent_second, + ] + _interfaces_form_alternating_quadrilateral(interfaces) || continue + allowed_choices = [ + [choice for (direction_index, choice) in records[label][interface] + if _joint_frame_choice_is_outward( + lattice, interface, directions[direction_index], interfaces, + )] + for (label, interface) in enumerate(interfaces) + ] + any(isempty, allowed_choices) && continue + push!(allowed_pairs, pair) + for choices in allowed_choices[2:4] + _sat_clause!(cnf, -opposite_variable, -pair, choices...) + end + end + _sat_clause!(cnf, -opposite_variable, allowed_pairs...) + end + + for choices in frame_choices, (vertex, direction_index, variable) in choices + pin = canonical[vertex] + direction = directions[direction_index] + for (site_vertex, site) in enumerate(canonical) + _joint_site_blocks_ray(site, pin, direction, directions) && + _sat_clause!(cnf, -variable, -selected[site_vertex]) + end + end + for first_label in 1:3, second_label in first_label+1:4, + (first_vertex, first_direction, first_variable) in frame_choices[first_label], + (second_vertex, second_direction, second_variable) in frame_choices[second_label] + first_pin = canonical[first_vertex] + second_pin = canonical[second_vertex] + first_ray = directions[first_direction] + second_ray = directions[second_direction] + first_start = (first_pin[1] + first_ray[1], first_pin[2] + first_ray[2]) + second_start = (second_pin[1] + second_ray[1], second_pin[2] + second_ray[2]) + _joint_rays_touch( + first_start, first_ray, second_start, second_ray, directions, + ) && _sat_clause!(cnf, -first_variable, -second_variable) + end + return frame_choices, choices_by_slot_vertex +end + +function _joint_crossing_sat_cnf( + target, lattice, shape, atom_count, offset; + canonical_shift=(0, 0), first_direction_index=1, +) + columns, rows = shape + coordinates = lattice isa Triangular ? sort!(_from_canonical.( + Ref(lattice), _LatticeCoordinate[ + (q + canonical_shift[1], r + canonical_shift[2]) + for q in -1:columns-2 for r in 1-rows:0 + ], + )) : _LatticeCoordinate[ + (column + canonical_shift[1], row + canonical_shift[2]) + for column in 0:columns-1 for row in 0:rows-1 + ] + canonical = _canonical_coordinate.(Ref(lattice), coordinates) + coordinate_index = Dict(point => vertex for (vertex, point) in enumerate(canonical)) + directions = _lattice_directions(lattice) + edge_list = Tuple{Int, Int}[] + for (vertex, point) in enumerate(canonical), direction in directions + other = get(coordinate_index, + (point[1] + direction[1], point[2] + direction[2]), 0) + vertex < other && push!(edge_list, (vertex, other)) + end + adjacent = [Int[] for _ in coordinates] + for (first, second) in edge_list + push!(adjacent[first], second) + push!(adjacent[second], first) + end + order, layers = _independent_set_bdd(coordinates, edge_list) + + cnf = _SatCnf() + selected = [_sat_variable!(cnf) for _ in coordinates] + _sat_exactly!(cnf, selected, atom_count) + origin = only(findall(==((0, 0)), canonical)) + _add_selected_connectivity!(cnf, selected, adjacent, origin, atom_count) + frame_choices, choices_by_slot_vertex = _add_joint_frame_choices!( + cnf, selected, coordinates, lattice, first_direction_index, + ) + for state in _essential_lower_states(target) + _add_joint_lower_state_constraint!( + cnf, selected, target, edge_list, frame_choices, state, offset, + ) + end + completion = _target_completion(target) + for state in _essential_upper_states(completion) + _add_joint_independent_set_upper_bound!( + cnf, selected, choices_by_slot_vertex, order, layers, state, + Int(completion[state+1] + offset), + ) + end + return cnf, selected, frame_choices, coordinates +end + +function _solve_joint_crossing_sat( + target, lattice, shape, atom_count, offset; + verbose=false, seconds=600, kissat_executable=nothing, initial_seed=1, + canonical_shift=(0, 0), first_direction_index=1, +) + cnf, selected, frame_choices, coordinates = _joint_crossing_sat_cnf( + target, lattice, shape, atom_count, offset; + canonical_shift, first_direction_index, + ) + solver = _new_sat_solver(cnf) + choice_variables = [choice[3] for choices in frame_choices for choice in choices] + projected = [selected; choice_variables] + model_index = 0 + deadline = time() + seconds + while time() < deadline + model_index += 1 + status, assignment = _next_joint_assignment!( + solver, projected; seed=initial_seed + model_index - 1, + seconds=max(1, floor(Int, deadline - time())), + kissat_executable, + ) + status == :unknown && return nothing + status == :unsat && return nothing + selected_assignment = assignment[1:length(selected)] + choice_assignment = assignment[length(selected)+1:end] + pins = _LatticeCoordinate[] + rays = Int[] + chosen_choices = Tuple{Int, Int, Int}[] + cursor = 0 + for choices in frame_choices + chosen = only(choice for choice in choices if choice_assignment[cursor += 1]) + push!(chosen_choices, chosen) + push!(pins, coordinates[chosen[1]]) + push!(rays, chosen[2]) + end + selected_indices = findall(identity, selected_assignment) + sites = coordinates[selected_indices] + checks = _check_crossing_frame(lattice, _LatticePatch(sites, pins, rays)) + if all(checks) + analysis = _analyze_crossing_candidate( + target, lattice, sites, pins, rays, + ) + analysis.solved && return analysis + error("joint SAT candidate failed a non-geometric encoded constraint") + end + chosen_ports = last.(chosen_choices) + intrinsic = _check_crossing_frame( + lattice, _LatticePatch(pins, pins, rays), + ) + clause = all(intrinsic) ? + [-selected[selected_indices]; -chosen_ports] : + -chosen_ports + _add_solver_clause!(solver, clause) + verbose && println((; model_index, frame_checks=checks, + intrinsic_frame_checks=intrinsic, blocker_length=length(clause))) + verbose && flush(stdout) + end + return nothing +end + +function _lower_state_filter_problem(target, context, atom_count, offset, state) + cnf = _SatCnf() + selected = [_sat_variable!(cnf) for _ in context.coordinates] + foreach(pin -> _sat_clause!(cnf, selected[pin]), context.boundary) + _sat_exactly!(cnf, selected, atom_count) + _add_lower_state_constraint!( + cnf, selected, target, context, state, offset, + ) + return _new_sat_solver(cnf), selected +end + +function _solve_next_fixed_crossing_sat!( + solver, selected, target, lattice, context, +) + assignment = _next_selected_assignment!(solver, selected) + assignment === nothing && return nothing + chosen = findall(identity, assignment) + sites = context.coordinates[chosen] + frame = context.frame + analysis = _analyze_crossing_candidate( + target, lattice, sites, frame.pins, frame.rays, + ) + analysis.solved || error("SAT candidate failed its encoded constraints") + return analysis +end + +mutable struct _KissatEnumerator + variables::Int + clauses::Vector{Vector{Int}} +end + +_new_sat_solver(cnf::_SatCnf) = + _KissatEnumerator(cnf.variables, copy(cnf.clauses)) +_add_solver_clause!(solver::_KissatEnumerator, clause) = + push!(solver.clauses, collect(clause)) +_kissat_init() = ccall((:kissat_init, Kissat_jll.libkissat), Ptr{Cvoid}, ()) +_kissat_add(solver, literal) = ccall((:kissat_add, Kissat_jll.libkissat), Cvoid, (Ptr{Cvoid}, Cint), solver, literal) +_kissat_solve(solver) = ccall((:kissat_solve, Kissat_jll.libkissat), Cint, (Ptr{Cvoid},), solver) +_kissat_value(solver, variable) = ccall((:kissat_value, Kissat_jll.libkissat), Cint, (Ptr{Cvoid}, Cint), solver, variable) +_kissat_release(solver) = ccall((:kissat_release, Kissat_jll.libkissat), Cvoid, (Ptr{Cvoid},), solver) +_kissat_quiet(solver) = ccall((:kissat_set_option, Kissat_jll.libkissat), Cint, (Ptr{Cvoid}, Cstring, Cint), solver, "quiet", 1) +_kissat_seed(solver, seed) = ccall( + (:kissat_set_option, Kissat_jll.libkissat), Cint, + (Ptr{Cvoid}, Cstring, Cint), solver, "seed", seed, +) +_kissat_walkinitially(solver) = ccall( + (:kissat_set_option, Kissat_jll.libkissat), Cint, + (Ptr{Cvoid}, Cstring, Cint), solver, "walkinitially", 1, +) +_kissat_set_conflict_limit(solver, limit) = ccall( + (:kissat_set_conflict_limit, Kissat_jll.libkissat), Cint, + (Ptr{Cvoid}, Cuint), solver, limit, +) + +function _next_selected_assignment!(enumerator::_KissatEnumerator, selected) + solver = _kissat_init() + _kissat_quiet(solver) + for clause in enumerator.clauses + foreach(literal -> _kissat_add(solver, literal), clause) + _kissat_add(solver, 0) + end + status = _kissat_solve(solver) + assignment = status == 10 ? + Bool[_kissat_value(solver, variable) > 0 for variable in selected] : nothing + _kissat_release(solver) + status == 20 && return nothing + status == 10 || error("SAT solver returned an undefined result") + push!(enumerator.clauses, [ + assignment[index] ? -variable : variable + for (index, variable) in enumerate(selected) + ]) + return assignment +end + +function _next_joint_assignment!( + enumerator::_KissatEnumerator, variables; + seed, seconds, kissat_executable, +) + path, stream = mktemp() + try + println(stream, "p cnf $(enumerator.variables) $(length(enumerator.clauses))") + for clause in enumerator.clauses + println(stream, join(clause, ' '), " 0") + end + close(stream) + output = IOBuffer() + executable = isnothing(kissat_executable) ? + Kissat_jll.kissat() : kissat_executable + command = `$executable --sat --walkinitially --seed=$seed -q --time=$seconds $path` + process = run(pipeline(ignorestatus(command), stdout=output, stderr=stderr)) + status = process.exitcode + status == 20 && return :unsat, nothing + status == 0 && return :unknown, nothing + status == 10 || error("Kissat exited with status $status") + positive = Set{Int}() + for line in eachline(seekstart(output)) + startswith(line, "v ") || continue + for literal in split(line)[2:end] + value = parse(Int, literal) + value > 0 && push!(positive, value) + end + end + return :sat, Bool[variable in positive for variable in variables] + finally + isopen(stream) && close(stream) + rm(path) + end +end + + +function _next_selected_assignment_limited!( + enumerator::_KissatEnumerator, selected, conflict_limit; seed=1, +) + solver = _kissat_init() + _kissat_quiet(solver) + _kissat_seed(solver, seed) + for clause in enumerator.clauses + foreach(literal -> _kissat_add(solver, literal), clause) + _kissat_add(solver, 0) + end + _kissat_set_conflict_limit(solver, conflict_limit) + status = _kissat_solve(solver) + assignment = status == 10 ? + Bool[_kissat_value(solver, variable) > 0 for variable in selected] : nothing + _kissat_release(solver) + status == 0 && return :unknown, nothing + status == 20 && return :unsat, nothing + status == 10 || error("SAT solver returned an undefined result") + push!(enumerator.clauses, [ + assignment[index] ? -variable : variable + for (index, variable) in enumerate(selected) + ]) + return :sat, assignment +end + +function _search_crossing_sat( + target_graph, target_boundary, target_reduced, lattice; + min_vertices, max_vertices, max_evaluations, max_frame_evaluations, + max_results, window_side, checkpoint_path, checkpoint_interval, +) + all(value -> isinf(value) || isinteger(value), target_reduced) || + error("the reduced alpha tensor must contain integers or -Inf") + completion = _target_completion(target_reduced) + lower_filter_states = _essential_lower_states(target_reduced) + resize!(lower_filter_states, min(2, length(lower_filter_states))) + gadgets = UnweightedGadget[] + seen_gadgets = Set{Tuple}() + stats = _SatSearchStats() + checkpoint = checkpoint_path === nothing || !isfile(checkpoint_path) ? + nothing : _read_sat_search_checkpoint(checkpoint_path) + if checkpoint !== nothing + checkpoint.target == Float64.(target_reduced) || + error("checkpoint target does not match this search") + checkpoint.lattice == _lattice_symbol(lattice) || + error("checkpoint lattice does not match this search") + checkpoint.min_vertices == min_vertices || + error("checkpoint min_vertices does not match this search") + checkpoint.max_vertices == max_vertices || + error("checkpoint max_vertices does not match this search") + checkpoint.window_side == window_side || + error("checkpoint window_side does not match this search") + stats = checkpoint.stats + end + evaluated = checkpoint === nothing ? 0 : checkpoint.evaluated + frame_evaluated = Ref(checkpoint === nothing ? 0 : checkpoint.frame_evaluated) + next_checkpoint_evaluation = evaluated + checkpoint_interval + for atom_count in min_vertices:max_vertices + checkpoint !== nothing && atom_count < checkpoint.atom_count && continue + atom_count > window_side^2 && break + seen_frames = Set{Tuple}() + offsets = _centered_offsets( + (-Int(minimum(completion))):(atom_count - Int(maximum(completion))), + ) + window_shapes = _crossing_window_shapes(window_side, atom_count) + for window_shape in window_shapes + if checkpoint !== nothing && atom_count == checkpoint.atom_count + shape_index = findfirst(==(checkpoint.window_shape), window_shapes) + shape_index === nothing && + error("checkpoint window shape is not part of this search") + current_index = findfirst(==(window_shape), window_shapes) + current_index < shape_index && continue + end + start_cursor = checkpoint !== nothing && + atom_count == checkpoint.atom_count && + window_shape == checkpoint.window_shape ? checkpoint.frame_cursor : 0 + start_offset = checkpoint !== nothing && + atom_count == checkpoint.atom_count && + window_shape == checkpoint.window_shape ? checkpoint.offset_cursor : 0 + start_order = checkpoint !== nothing && + atom_count == checkpoint.atom_count && + window_shape == checkpoint.window_shape ? checkpoint.order_cursor : 0 + save_checkpoint = (cursor, order, offset) -> begin + checkpoint_path === nothing && return + evaluated >= next_checkpoint_evaluation || return + _write_sat_search_checkpoint( + checkpoint_path, + _SatSearchCheckpoint( + Float64.(target_reduced), _lattice_symbol(lattice), + min_vertices, max_vertices, window_side, atom_count, + window_shape, cursor, order, offset, evaluated, + frame_evaluated[], stats, + ), + ) + next_checkpoint_evaluation = evaluated + checkpoint_interval + end + stopped = _foreach_crossing_frame_clique( + lattice, window_shape, frame_evaluated, max_frame_evaluations; + min_allowed=atom_count, seen=seen_frames, start_cursor, + start_order, + ) do frame, frame_cursor, order_cursor + stats.frame_candidates += 1 + isempty(offsets) && return nothing + evaluated == max_evaluations && return :budget + context = _prepare_sat_frame(lattice, frame) + resume_offset = frame_cursor == start_cursor && + order_cursor == start_order ? start_offset : 0 + for (offset_index, offset) in enumerate(offsets) + offset_index <= resume_offset && continue + save_checkpoint(frame_cursor, order_cursor, offset_index - 1) + rejected = false + for (filter_index, state) in enumerate(lower_filter_states) + evaluated == max_evaluations && return :budget + filter_solver, filter_selected = _lower_state_filter_problem( + target_reduced, context, atom_count, offset, state, + ) + evaluated += 1 + save_checkpoint(frame_cursor, order_cursor, offset_index - 1) + if _next_selected_assignment!( + filter_solver, filter_selected, + ) === nothing + filter_index == 1 ? + (stats.first_lower_rejected += 1) : + (stats.second_lower_rejected += 1) + rejected = true + break + end + end + if rejected + save_checkpoint(frame_cursor, order_cursor, offset_index) + continue + end + evaluated == max_evaluations && return :budget + solver, selected = _fixed_crossing_sat_problem( + target_reduced, context, atom_count, offset, + ) + stats.full_solves += 1 + while true + evaluated == max_evaluations && return :budget + evaluated += 1 + save_checkpoint(frame_cursor, order_cursor, offset_index - 1) + analysis = _solve_next_fixed_crossing_sat!( + solver, selected, target_reduced, lattice, context, + ) + analysis === nothing && break + valid, verified_offset = is_gadget_replacement( + target_graph, analysis.graph, target_boundary, + analysis.boundary, + ) + valid || error("SAT candidate failed the fixed verifier") + analysis.offset == verified_offset || + error("SAT candidate offset mismatch") + key = _canonical_crossing_frame_key(lattice, ( + allowed=analysis.patch.coordinates, + pins=analysis.patch.pins, + rays=analysis.patch.rays, + )) + key in seen_gadgets && continue + push!(seen_gadgets, key) + push!( + gadgets, + _unweighted_gadget(target_graph, lattice, analysis), + ) + length(gadgets) == max_results && return :solution + end + save_checkpoint(frame_cursor, order_cursor, offset_index) + end + return nothing + end + stopped === nothing || return gadgets, evaluated, stopped + end + end + return gadgets, evaluated, :search_space_exhausted +end + +function _crossing_window_shapes(side, atom_count) + minimum_minor = cld(atom_count, side) + shapes = [(side, side)] + for minor in side-1:-1:minimum_minor + push!(shapes, (side, minor)) + push!(shapes, (minor, side)) + end + return shapes +end + +function _single_pin_corridor_clear(lattice, pins, index, ray_index) + direction = _lattice_directions(lattice)[ray_index] + interface = _lattice_step(lattice, pins[index], direction, 1) + start = _canonical_coordinate(lattice, interface) + canonical_pins = _canonical_coordinate.(Ref(lattice), pins) + adjacency_offsets = [_LatticeCoordinate[(0, 0)]; _lattice_directions(lattice)] + for pin_index in eachindex(pins) + pin_index == index && continue + site = canonical_pins[pin_index] + for offset in adjacency_offsets + _point_on_ray( + (site[1] + offset[1], site[2] + offset[2]), start, direction, + ) && return false + end + end + return true +end + +function _pin_rays_are_compatible(lattice, pins, rays) + directions = _lattice_directions(lattice)[collect(rays)] + interfaces = [ + _lattice_step(lattice, pin, direction, 1) + for (pin, direction) in zip(pins, directions) + ] + starts = _canonical_coordinate.(Ref(lattice), interfaces) + adjacency_offsets = [_LatticeCoordinate[(0, 0)]; _lattice_directions(lattice)] + for first in 1:3, second in first+1:4, offset in adjacency_offsets + _rays_touch( + starts[first], directions[first], starts[second], directions[second], + offset, + ) && return false + end + return true +end + +function _labels_alternate(lattice, pins, rays) + directions = _lattice_directions(lattice)[rays] + interfaces = [ + _lattice_step(lattice, pin, direction, 1) + for (pin, direction) in zip(pins, directions) + ] + interface_geometry = _geometry_coordinate.(Ref(lattice), interfaces) + occupied_geometry = _geometry_coordinate.(Ref(lattice), pins) + hull = _strict_convex_hull([occupied_geometry; interface_geometry]) + return _interfaces_alternate(hull, interface_geometry) +end + +function _allowed_crossing_sites(lattice, window, pins, rays, min_allowed) + directions = _lattice_directions(lattice)[rays] + interfaces = [ + _lattice_step(lattice, pin, direction, 1) + for (pin, direction) in zip(pins, directions) + ] + interface_geometry = _geometry_coordinate.(Ref(lattice), interfaces) + pin_geometry = _geometry_coordinate.(Ref(lattice), pins) + starts = _canonical_coordinate.(Ref(lattice), interfaces) + adjacency_offsets = [_LatticeCoordinate[(0, 0)]; _lattice_directions(lattice)] + allowed = copy(pins) + candidates = [site for site in window if site ∉ pins] + for (candidate_index, site) in enumerate(candidates) + canonical_site = _canonical_coordinate(lattice, site) + corridor_clear = all(eachindex(starts)) do ray_index + all(adjacency_offsets) do offset + !_point_on_ray( + (canonical_site[1] + offset[1], canonical_site[2] + offset[2]), + starts[ray_index], directions[ray_index], + ) + end + end + if corridor_clear + site_geometry = _geometry_coordinate(lattice, site) + hull = _strict_convex_hull([ + pin_geometry; site_geometry; interface_geometry + ]) + all(in(hull), interface_geometry) && push!(allowed, site) + end + remaining = length(candidates) - candidate_index + length(allowed) + remaining >= min_allowed || return nothing + end + sort!(allowed) + return allowed +end + +function _crossing_port_candidates(lattice, window) + return [ + (pin, ray) for pin in window + for ray in eachindex(_lattice_directions(lattice)) + ] +end + +function _crossing_ports_compatible(lattice, first, second) + first_pin, first_ray = first + second_pin, second_ray = second + first_pin == second_pin && return false + pins = [first_pin, second_pin] + _single_pin_corridor_clear(lattice, pins, 1, first_ray) || return false + _single_pin_corridor_clear(lattice, pins, 2, second_ray) || return false + directions = _lattice_directions(lattice) + first_direction = directions[first_ray] + second_direction = directions[second_ray] + first_start = _canonical_coordinate( + lattice, _lattice_step(lattice, first_pin, first_direction, 1), + ) + second_start = _canonical_coordinate( + lattice, _lattice_step(lattice, second_pin, second_direction, 1), + ) + return all([_LatticeCoordinate[(0, 0)]; directions]) do offset + !_rays_touch( + first_start, first_direction, second_start, second_direction, offset, + ) + end +end + +function _crossing_port_prefix(lattice, window, rank) + candidates = _crossing_port_candidates(lattice, window) + compatible_rank = 0 + for first in 1:length(candidates)-1, second in first+1:length(candidates) + _crossing_ports_compatible( + lattice, candidates[first], candidates[second], + ) || continue + compatible_rank += 1 + compatible_rank == rank && return (first, second) + end + error("compatible port-prefix rank $rank exceeds $compatible_rank") +end + +function _foreach_crossing_port_clique( + visit, lattice, window; shard_index=0, shard_count=1, + port_prefix=nothing, +) + candidates = _crossing_port_candidates(lattice, window) + compatible = falses(length(candidates), length(candidates)) + for first in 1:length(candidates)-1, second in first+1:length(candidates) + compatible[first, second] = _crossing_ports_compatible( + lattice, candidates[first], candidates[second], + ) + end + chosen = Int[] + function extend(available) + if length(chosen) == 4 + return visit(candidates[chosen]) + end + needed = 4 - length(chosen) + for position in 1:length(available)-needed+1 + candidate = available[position] + if isempty(chosen) + port_prefix !== nothing && candidate != port_prefix[1] && continue + elseif length(chosen) == 1 + port_prefix !== nothing && candidate != port_prefix[2] && continue + pair_key = (chosen[1] - 1) * length(candidates) + candidate - 1 + mod(pair_key, shard_count) == shard_index || continue + end + push!(chosen, candidate) + remaining = [ + other for other in @view(available[position+1:end]) + if compatible[candidate, other] + ] + result = extend(remaining) + pop!(chosen) + result === nothing || return result + end + return nothing + end + return extend(collect(eachindex(candidates))) +end + +function _foreach_crossing_frame_clique( + visit, lattice, shape, evaluated, limit; + min_allowed=0, seen=Set{Tuple}(), start_cursor=0, start_order=0, + shard_index=0, shard_count=1, port_prefix=nothing, +) + columns, rows = shape + window = _LatticeCoordinate[ + (column, row) for column in 0:columns-1 for row in 0:rows-1 + ] + cursor = 0 + return _foreach_crossing_port_clique( + lattice, window; shard_index, shard_count, port_prefix, + ) do ports + cursor += 1 + cursor < start_cursor && return nothing + evaluated[] == limit && return :frame_budget + evaluated[] += 1 + pins = first.(ports) + rays = last.(ports) + checks = _check_crossing_frame( + lattice, _LatticePatch(pins, pins, rays), + ) + checks[1] && checks[3] && checks[4] || return nothing + directions = _lattice_directions(lattice)[rays] + interfaces = [ + _lattice_step(lattice, pin, direction, 1) + for (pin, direction) in zip(pins, directions) + ] + interface_geometry = _geometry_coordinate.(Ref(lattice), interfaces) + occupied_geometry = _geometry_coordinate.(Ref(lattice), pins) + hull = _strict_convex_hull([occupied_geometry; interface_geometry]) + valid_orders = [ + collect(order) for order in permutations(1:4) + if _interfaces_alternate(hull, interface_geometry[collect(order)]) + ] + isempty(valid_orders) && return nothing + allowed = _allowed_crossing_sites( + lattice, window, pins, rays, min_allowed, + ) + allowed === nothing && return nothing + for (order_cursor, order) in enumerate(valid_orders) + cursor == start_cursor && order_cursor < start_order && continue + frame = ( + pins=pins[order], rays=rays[order], allowed=allowed, + ) + key = _canonical_crossing_frame_key(lattice, frame) + key in seen && continue + push!(seen, key) + result = visit(frame, cursor, order_cursor) + result === nothing || return result + end + return nothing + end +end diff --git a/src/core/unweighted_search.jl b/src/core/unweighted_search.jl index dd28596..d74004b 100644 --- a/src/core/unweighted_search.jl +++ b/src/core/unweighted_search.jl @@ -1,90 +1,800 @@ # ============================================================================ -# Unweighted Gadget Types +# Unweighted lattice search # ============================================================================ -""" - UnweightedGadget +const _LatticeCoordinate = Tuple{Int, Int} -Result of an unweighted gadget search. -Stores the pattern graph R, replacement graph R', boundary vertices, -constant offset between reduced alpha tensors, and optional vertex positions. -""" +"""A verifier-accepted replacement together with its concrete lattice embedding.""" struct UnweightedGadget pattern_graph::SimpleGraph{Int} replacement_graph::SimpleGraph{Int} boundary_vertices::Vector{Int} constant_offset::Float64 - pos::Union{Nothing, Vector{Tuple{Float64, Float64}}} + lattice::Symbol + lattice_coordinates::Vector{_LatticeCoordinate} + pos::Vector{Tuple{Float64, Float64}} + pin_rays::Vector{_LatticeCoordinate} end -# ============================================================================ -# Unweighted Filter Construction -# ============================================================================ +"""Outcome of a bounded direct SAT search on one concrete lattice.""" +struct UnweightedSearchResult + target_graph::SimpleGraph{Int} + target_boundary::Vector{Int} + lattice::Symbol + gadgets::Vector{UnweightedGadget} + evaluated::Int + termination_reason::Symbol +end + +"""One verifier-certified atom-reducing rewrite.""" +struct UnweightedRewriteStep + rule::Symbol + before::UnweightedGadget + after::UnweightedGadget +end + +"""Result of rewriting and re-synthesizing one realized lattice gadget.""" +struct UnweightedOptimizationResult + gadget::UnweightedGadget + steps::Vector{UnweightedRewriteStep} + sat_evaluations::Int + unresolved_sat_evaluations::Int + termination_reason::Symbol +end +struct _LatticePatch + coordinates::Vector{_LatticeCoordinate} + pins::Vector{_LatticeCoordinate} + rays::Vector{Int} +end + +""" +Search a four-pin unweighted gadget in a finite window of `lattice`. +Every returned gadget satisfies the reduced-alpha target up to a constant and +the four-direction crossing geometry. """ - _make_unweighted_filter(pattern_graph, pattern_boundary; prefilter) +function search_unweighted_gadgets( + target_graph::SimpleGraph{Int}, + target_boundary::Vector{Int}, + lattice::LatticeType=Triangular(); + min_vertices::Int=length(target_boundary) + 1, + max_vertices::Int=min_vertices + 8, + max_evaluations::Int=2_000, + max_frame_evaluations::Int=1_000_000, + max_results::Int=1, + window_side::Int=4, + checkpoint_path::Union{Nothing, String}=nothing, + checkpoint_interval::Int=10_000, +) + boundary_count = length(target_boundary) + boundary_count == 4 || + throw(ArgumentError("unweighted lattice search requires four boundary vertices")) + min_vertices >= boundary_count || + throw(ArgumentError("min_vertices must be at least the number of boundary vertices")) + max_vertices >= min_vertices || + throw(ArgumentError("max_vertices must be at least min_vertices")) + max_evaluations > 0 || + throw(ArgumentError("max_evaluations must be positive")) + max_frame_evaluations > 0 || + throw(ArgumentError("max_frame_evaluations must be positive")) + max_results > 0 || throw(ArgumentError("max_results must be positive")) + window_side >= 2 || throw(ArgumentError("window_side must be at least 2")) + checkpoint_interval > 0 || + throw(ArgumentError("checkpoint_interval must be positive")) + + target_reduced = vec(calculate_reduced_alpha_tensor(target_graph, target_boundary)) + all(isinf, target_reduced) && + error("target graph has an entirely -Inf reduced alpha tensor") + gadgets, evaluated, reason = _search_crossing_sat( + target_graph, target_boundary, target_reduced, lattice; + min_vertices, max_vertices, max_evaluations, max_frame_evaluations, + max_results, window_side, checkpoint_path, checkpoint_interval, + ) + return UnweightedSearchResult( + target_graph, copy(target_boundary), _lattice_symbol(lattice), + gadgets, evaluated, reason, + ) +end -Build a filter closure that checks candidate graphs against the target pattern. """ -function _make_unweighted_filter( - pattern_graph::SimpleGraph{Int}, - pattern_boundary::Vector{Int}; - prefilter::Bool=true, +Search one joint occupancy-and-frame SAT instance in a fixed lattice window. + +Unlike `search_unweighted_gadgets`, this formulation chooses the occupied sites, +four pins, and four outward rays in one CNF. `atom_count` and `offset` identify +the exact instance to solve. The unchanged gadget verifier checks every result. +""" +function search_unweighted_gadget_joint( + target_graph::SimpleGraph{Int}, + target_boundary::Vector{Int}, + lattice::LatticeType=Triangular(); + window_shape::Tuple{Int, Int}, + atom_count::Int, + offset::Int, + seconds::Int=600, + kissat_executable::Union{Nothing, String}=nothing, + seed::Int=1, + canonical_shift::Tuple{Int, Int}=(0, 0), + first_ray::Int=1, ) - k = length(pattern_boundary) - target_reduced = vec(calculate_reduced_alpha_tensor(pattern_graph, pattern_boundary)) - all(isinf, target_reduced) && error("target graph has an entirely -Inf reduced alpha tensor") - target_mask = inf_mask(target_reduced) - apply_prefilter = prefilter && pins_prefilter(pattern_graph, pattern_boundary) - return function(candidate::SimpleGraph{Int}, pos, pin_set) - vertex_pool = something(pin_set, 1:Graphs.nv(candidate)) - if apply_prefilter && !pins_prefilter(candidate, vertex_pool) - return nothing + length(target_boundary) == 4 || + throw(ArgumentError("joint unweighted search requires four boundary vertices")) + all(>=(2), window_shape) || + throw(ArgumentError("window dimensions must be at least 2")) + length(target_boundary) <= atom_count <= prod(window_shape) || + throw(ArgumentError("atom_count must fit in the window")) + seconds > 0 || throw(ArgumentError("seconds must be positive")) + seed > 0 || throw(ArgumentError("seed must be positive")) + directions = _lattice_directions(lattice) + first_ray in eachindex(directions) || + throw(ArgumentError("first_ray is not a lattice direction index")) + + target_reduced = vec(calculate_reduced_alpha_tensor( + target_graph, target_boundary, + )) + all(value -> isinf(value) || isinteger(value), target_reduced) || + error("the reduced alpha tensor must contain integers or -Inf") + analysis = _solve_joint_crossing_sat( + target_reduced, lattice, window_shape, atom_count, offset; + seconds, kissat_executable, initial_seed=seed, canonical_shift, + first_direction_index=first_ray, + ) + analysis === nothing && return nothing + verified, verified_offset = is_gadget_replacement( + target_graph, analysis.graph, target_boundary, analysis.boundary, + ) + verified || error("joint SAT result failed the gadget verifier") + verified_offset == analysis.offset || + error("joint SAT result and gadget verifier disagree on the offset") + return _unweighted_gadget(target_graph, lattice, analysis) +end + +""" + optimize_unweighted_gadget(gadget, target_boundary; kwargs...) + +Reduce a verifier-accepted four-pin lattice gadget with certified rewrite rules. +The optimizer contracts even boundary tails, contracts opposite leaf pins in +pairs, and uses fixed-frame SAT to re-synthesize the interior after a one-step +frame rewrite. Direct rewrites are explored to a closure before re-synthesis, +so the optimizer does not commit to the first smaller direct result. It never +proposes arbitrary vertex deletion. +""" +function optimize_unweighted_gadget( + gadget::UnweightedGadget, + target_boundary::Vector{Int}; + min_vertices::Int=length(target_boundary), + max_sat_evaluations::Int=256, + max_sat_conflicts::Int=100_000, + host_radius::Int=1, +) + length(target_boundary) == 4 || + throw(ArgumentError("unweighted rewrite optimization requires four target boundary vertices")) + min_vertices >= length(target_boundary) || + throw(ArgumentError("min_vertices must include all boundary vertices")) + min_vertices <= nv(gadget.replacement_graph) || + throw(ArgumentError("min_vertices exceeds the current gadget size")) + max_sat_evaluations >= 0 || + throw(ArgumentError("max_sat_evaluations must be nonnegative")) + max_sat_conflicts >= 0 || + throw(ArgumentError("max_sat_conflicts must be nonnegative")) + host_radius >= 0 || throw(ArgumentError("host_radius must be nonnegative")) + + lattice = _gadget_lattice(gadget) + target_reduced = vec(calculate_reduced_alpha_tensor( + gadget.pattern_graph, target_boundary, + )) + valid, offset = is_gadget_replacement( + gadget.pattern_graph, gadget.replacement_graph, target_boundary, + gadget.boundary_vertices, + ) + valid || throw(ArgumentError("the input gadget does not replace the requested target boundary")) + offset == gadget.constant_offset || + throw(ArgumentError("the input gadget stores the wrong constant offset")) + current = gadget + steps = UnweightedRewriteStep[] + sat_evaluations = 0 + unresolved_sat_evaluations = 0 + while true + nv(current.replacement_graph) == min_vertices && + return UnweightedOptimizationResult( + current, steps, sat_evaluations, unresolved_sat_evaluations, + :minimum_vertices, + ) + + current_key = _rewrite_state_key(current) + direct_states = [(gadget=current, steps=steps)] + seen_states = Set{Tuple}([current_key]) + for state in direct_states + for candidate in _direct_unweighted_rewrites( + state.gadget, target_boundary, target_reduced, lattice, min_vertices, + ) + candidate_steps = [ + state.steps; + UnweightedRewriteStep( + candidate.rule, state.gadget, candidate.gadget, + ) + ] + nv(candidate.gadget.replacement_graph) == min_vertices && + return UnweightedOptimizationResult( + candidate.gadget, candidate_steps, sat_evaluations, + unresolved_sat_evaluations, + :minimum_vertices, + ) + candidate_key = _rewrite_state_key(candidate.gadget) + candidate_key in seen_states && continue + push!(seen_states, candidate_key) + push!(direct_states, (; + gadget=candidate.gadget, steps=candidate_steps, + )) + end end - (Graphs.nv(candidate) < k || length(vertex_pool) < k) && return nothing - for boundary in Combinatorics.combinations(vertex_pool, k) - candidate_reduced = vec(calculate_reduced_alpha_tensor(candidate, boundary)) - all(isinf, candidate_reduced) && continue - candidate_mask = inf_mask(candidate_reduced) - candidate_mask == target_mask || continue - valid, constant_offset = is_diff_by_constant(candidate_reduced, target_reduced) - if valid - return UnweightedGadget( - pattern_graph, candidate, boundary, - constant_offset, pos) + + best = argmin( + state -> nv(state.gadget.replacement_graph), direct_states, + ) + sources = [ + sort(direct_states[2:end]; by=state -> -nv(state.gadget.replacement_graph)); + direct_states[1] + ] + resumed = false + for source in sources + candidate, used, unresolved, status = _resynthesized_unweighted_rewrite( + source.gadget, target_boundary, target_reduced, lattice, min_vertices, + max_sat_evaluations - sat_evaluations, max_sat_conflicts, + host_radius, + ) + sat_evaluations += used + unresolved_sat_evaluations += unresolved + if candidate !== nothing + current = candidate.gadget + steps = [ + source.steps; + UnweightedRewriteStep( + candidate.rule, source.gadget, candidate.gadget, + ) + ] + resumed = true + break end + status == :budget && return UnweightedOptimizationResult( + best.gadget, best.steps, sat_evaluations, + unresolved_sat_evaluations, :sat_budget, + ) end - return nothing + resumed || return UnweightedOptimizationResult( + best.gadget, best.steps, sat_evaluations, + unresolved_sat_evaluations, + iszero(unresolved_sat_evaluations) ? + :rewrite_fixed_point : :sat_unknown, + ) end end -# ============================================================================ -# Unweighted Search -# ============================================================================ +_rewrite_state_key(gadget) = ( + Tuple(gadget.lattice_coordinates), + Tuple(gadget.lattice_coordinates[gadget.boundary_vertices]), + Tuple(gadget.pin_rays), +) + +_gadget_lattice(gadget::UnweightedGadget) = + gadget.lattice == :KSG ? Square() : + gadget.lattice == :triangular ? Triangular() : + error("unknown gadget lattice $(gadget.lattice)") + +function _gadget_patch(gadget, lattice) + pins = gadget.lattice_coordinates[gadget.boundary_vertices] + directions = _lattice_directions(lattice) + rays = [_lattice_direction_index(directions, ray) for ray in gadget.pin_rays] + return _LatticePatch(copy(gadget.lattice_coordinates), pins, rays) +end + +function _rewrite_gadget( + gadget, target_boundary, target_reduced, lattice, sites, pins, rays, +) + analysis = _analyze_crossing_candidate( + target_reduced, lattice, collect(sites), collect(pins), collect(rays), + ) + analysis.solved || return nothing + return _certified_rewrite_gadget(gadget, target_boundary, lattice, analysis) +end + +function _certified_rewrite_gadget(gadget, target_boundary, lattice, analysis) + valid, offset = is_gadget_replacement( + gadget.pattern_graph, analysis.graph, target_boundary, analysis.boundary, + ) + valid || return nothing + offset == analysis.offset || error("rewrite analysis and verifier disagree on the offset") + return _unweighted_gadget(gadget.pattern_graph, lattice, analysis) +end + +function _direct_unweighted_rewrites( + gadget, target_boundary, target_reduced, lattice, min_vertices, +) + patch = _gadget_patch(gadget, lattice) + graph = gadget.replacement_graph + boundary = gadget.boundary_vertices + rewrites = NamedTuple[] + + nv(graph) - 2 >= min_vertices && for label in eachindex(boundary) + pin = boundary[label] + degree(graph, pin) == 1 || continue + middle = only(neighbors(graph, pin)) + for endpoint in neighbors(graph, middle) + endpoint == pin && continue + endpoint in boundary && continue + sites = setdiff(patch.coordinates, patch.coordinates[[pin, middle]]) + pins = copy(patch.pins) + pins[label] = patch.coordinates[endpoint] + rays = copy(patch.rays) + rays[label] = _direction_index_between( + lattice, pins[label], patch.coordinates[middle], + ) + rewritten = _rewrite_gadget( + gadget, target_boundary, target_reduced, lattice, sites, pins, rays, + ) + rewritten === nothing || push!(rewrites, (; + rule=:even_boundary_tail_contraction, gadget=rewritten, + )) + end + end + + nv(graph) - 2 >= min_vertices && for (first_label, second_label) in ((1, 3), (2, 4)) + first_pin = boundary[first_label] + second_pin = boundary[second_label] + degree(graph, first_pin) == 1 || continue + degree(graph, second_pin) == 1 || continue + first_neighbor = only(neighbors(graph, first_pin)) + second_neighbor = only(neighbors(graph, second_pin)) + first_neighbor == second_neighbor && continue + first_neighbor in boundary && continue + second_neighbor in boundary && continue + sites = setdiff( + patch.coordinates, patch.coordinates[[first_pin, second_pin]], + ) + pins = copy(patch.pins) + pins[first_label] = patch.coordinates[first_neighbor] + pins[second_label] = patch.coordinates[second_neighbor] + rays = copy(patch.rays) + rays[first_label] = _direction_index_between( + lattice, pins[first_label], patch.coordinates[first_pin], + ) + rays[second_label] = _direction_index_between( + lattice, pins[second_label], patch.coordinates[second_pin], + ) + rewritten = _rewrite_gadget( + gadget, target_boundary, target_reduced, lattice, sites, pins, rays, + ) + rewritten === nothing || push!(rewrites, (; + rule=:opposite_leaf_pin_contraction, gadget=rewritten, + )) + end + return rewrites +end + +function _direction_index_between(lattice, source, destination) + source_canonical = _canonical_coordinate(lattice, source) + destination_canonical = _canonical_coordinate(lattice, destination) + direction = ( + destination_canonical[1] - source_canonical[1], + destination_canonical[2] - source_canonical[2], + ) + return _lattice_direction_index(_lattice_directions(lattice), direction) +end + +function _expanded_lattice_host(lattice, coordinates, radius) + host = Set(coordinates) + frontier = Set(coordinates) + directions = _lattice_directions(lattice) + for _ in 1:radius + next_frontier = Set{_LatticeCoordinate}() + for coordinate in frontier, direction in directions + neighbor = _lattice_step(lattice, coordinate, direction, 1) + neighbor in host || push!(next_frontier, neighbor) + end + union!(host, next_frontier) + frontier = next_frontier + end + return sort!(collect(host)) +end + +function _one_step_frame_rewrites(lattice, patch) + directions = _lattice_directions(lattice) + frames = NamedTuple[] + for label in eachindex(patch.pins), step in directions + pins = copy(patch.pins) + pins[label] = _lattice_step(lattice, pins[label], step, 1) + length(unique(pins)) == 4 || continue + for ray in eachindex(directions) + rays = copy(patch.rays) + rays[label] = ray + _pin_rays_are_compatible(lattice, pins, rays) || continue + _labels_alternate(lattice, pins, rays) || continue + push!(frames, (; label, pins, rays)) + end + end + return frames +end + +function _resynthesized_unweighted_rewrite( + gadget, target_boundary, target_reduced, lattice, min_vertices, + sat_budget, max_sat_conflicts, host_radius, +) + iszero(sat_budget) && return nothing, 0, 0, :budget + patch = _gadget_patch(gadget, lattice) + window = _expanded_lattice_host(lattice, patch.coordinates, host_radius) + rewritten_frames = _one_step_frame_rewrites(lattice, patch) + sort!(rewritten_frames; by=frame -> ( + degree(gadget.replacement_graph, gadget.boundary_vertices[frame.label]) == 1, + frame.rays[frame.label] != patch.rays[frame.label], + frame.label, + )) + frame_hosts = [begin + allowed = _allowed_crossing_sites( + lattice, window, frame.pins, frame.rays, min_vertices, + ) + allowed === nothing ? nothing : (; frame, allowed) + end for frame in rewritten_frames] + filter!(!isnothing, frame_hosts) + contexts = Dict{Int, _SatFrameContext}() + completion = _target_completion(target_reduced) + evaluated = 0 + unresolved = 0 + for atom_count in nv(gadget.replacement_graph)-1:-1:min_vertices + expected_offset = Int(gadget.constant_offset) - + (nv(gadget.replacement_graph) - atom_count) + offsets = _centered_offsets( + -Int(minimum(completion)): + atom_count-Int(maximum(completion)), expected_offset, + ) + for offset in offsets + for (frame_index, frame_host) in enumerate(frame_hosts) + evaluated == sat_budget && + return nothing, evaluated, unresolved, :budget + length(frame_host.allowed) >= atom_count || continue + frame = (; + pins=frame_host.frame.pins, + rays=frame_host.frame.rays, + allowed=frame_host.allowed, + ) + context = get!(contexts, frame_index) do + _prepare_sat_frame(lattice, frame) + end + evaluated += 1 + solver, selected = _fixed_crossing_sat_problem( + target_reduced, context, atom_count, offset, + ) + status, assignment = _next_selected_assignment_limited!( + solver, selected, max_sat_conflicts; seed=evaluated, + ) + if status == :unknown + unresolved += 1 + continue + end + status == :unsat && continue + chosen = findall(identity, assignment) + sites = context.coordinates[chosen] + analysis = _analyze_crossing_candidate( + target_reduced, lattice, sites, frame.pins, frame.rays, + ) + analysis.solved || + error("SAT candidate failed its encoded constraints") + rewritten = _certified_rewrite_gadget( + gadget, target_boundary, lattice, analysis, + ) + rewritten === nothing && + error("fixed-frame SAT result failed the gadget verifier") + return (; + rule=:frame_rewrite_resynthesis, + gadget=rewritten, + ), evaluated, unresolved, :found + end + end + end + return nothing, evaluated, unresolved, :exhausted +end + +function _analyze_crossing_candidate(target_reduced, lattice, sites, pins, rays) + patch = _LatticePatch(sort(copy(sites)), copy(pins), copy(rays)) + graph, boundary, positions = _materialize_lattice_patch(lattice, patch) + reduced = vec(calculate_reduced_alpha_tensor(graph, boundary)) + tensor_valid, offset = is_diff_by_constant(reduced, target_reduced) + geometry_valid = all(_check_crossing_frame(lattice, patch)) + solved = tensor_valid && geometry_valid && is_connected(graph) + return (; + graph, boundary, positions, patch, solved, offset=Float64(offset), + ) +end + +function _unweighted_gadget(target_graph, lattice, analysis) + return UnweightedGadget( + target_graph, analysis.graph, analysis.boundary, analysis.offset, + _lattice_symbol(lattice), copy(analysis.patch.coordinates), + analysis.positions, _patch_ray_directions(lattice, analysis.patch), + ) +end + +_rotate_lattice_coordinate(::Triangular, point) = (-point[2], point[1] + point[2]) +_reflect_lattice_coordinate(::Triangular, point) = (point[2], point[1]) +_rotate_lattice_coordinate(::Square, point) = (-point[2], point[1]) +_reflect_lattice_coordinate(::Square, point) = (point[1], -point[2]) + +function _transform_lattice_coordinate(lattice, point, rotations, reflected) + transformed = reflected ? _reflect_lattice_coordinate(lattice, point) : point + for _ in 1:rotations + transformed = _rotate_lattice_coordinate(lattice, transformed) + end + return transformed +end + +function _canonical_crossing_frame_key(lattice, frame) + allowed = _canonical_coordinate.(Ref(lattice), frame.allowed) + pins = _canonical_coordinate.(Ref(lattice), frame.pins) + rays = _lattice_directions(lattice)[frame.rays] + return minimum((begin + transformed_allowed = _transform_lattice_coordinate.( + Ref(lattice), allowed, rotations, reflected, + ) + transformed_pins = _transform_lattice_coordinate.( + Ref(lattice), pins, rotations, reflected, + ) + transformed_rays = _transform_lattice_coordinate.( + Ref(lattice), rays, rotations, reflected, + ) + minimum_first = minimum(first, transformed_allowed) + minimum_last = minimum(last, transformed_allowed) + normalize(point) = (point[1] - minimum_first, point[2] - minimum_last) + ( + Tuple(sort(normalize.(transformed_allowed))), + Tuple(normalize.(transformed_pins)), + Tuple(transformed_rays), + ) + end for reflected in (false, true) for rotations in + 0:(lattice isa Triangular ? 5 : 3))) +end + +function _materialize_lattice_patch(lattice::LatticeType, patch::_LatticePatch) + positions = get_physical_positions(lattice, patch.coordinates) + graph = unit_disk_graph(positions, get_radius(lattice)) + coordinate_index = Dict(coordinate => index for (index, coordinate) in enumerate(patch.coordinates)) + boundary = [coordinate_index[pin] for pin in patch.pins] + return graph, boundary, positions +end + +function _normalize_lattice_patch( + ::Square, + coordinates::Vector{_LatticeCoordinate}, + pins::Vector{_LatticeCoordinate}, + rays::Vector{Int}, +) + min_x = minimum(first, coordinates) + min_y = minimum(last, coordinates) + translate(point) = (point[1] - min_x, point[2] - min_y) + return _LatticePatch(sort!(translate.(coordinates)), translate.(pins), copy(rays)) +end + +function _normalize_lattice_patch( + ::Triangular, + coordinates::Vector{_LatticeCoordinate}, + pins::Vector{_LatticeCoordinate}, + rays::Vector{Int}, +) + axial = [_offset_to_axial(point) for point in coordinates] + pin_axial = [_offset_to_axial(point) for point in pins] + min_q = minimum(first, axial) + min_r = minimum(last, axial) + translate(point) = (point[1] - min_q, point[2] - min_r) + translated = _axial_to_offset.(translate.(axial)) + translated_pins = _axial_to_offset.(translate.(pin_axial)) + return _LatticePatch(sort!(translated), translated_pins, copy(rays)) +end + +_offset_to_axial(point::_LatticeCoordinate) = (point[1] - fld(point[2], 2), point[2]) +_axial_to_offset(point::_LatticeCoordinate) = (point[1] + fld(point[2], 2), point[2]) + +_lattice_directions(::Square) = _LatticeCoordinate[ + (-1, -1), (0, -1), (1, -1), (-1, 0), (1, 0), (-1, 1), (0, 1), (1, 1), +] +_lattice_directions(::Triangular) = _LatticeCoordinate[ + (1, 0), (0, 1), (-1, 1), (-1, 0), (0, -1), (1, -1), +] + +function _lattice_step(::Square, point::_LatticeCoordinate, direction::_LatticeCoordinate, distance::Int) + return (point[1] + distance * direction[1], point[2] + distance * direction[2]) +end + +function _lattice_step(::Triangular, point::_LatticeCoordinate, direction::_LatticeCoordinate, distance::Int) + q, r = _offset_to_axial(point) + return _axial_to_offset((q + distance * direction[1], r + distance * direction[2])) +end + +function _lattice_symbol(lattice::LatticeType) + return lattice isa Square ? :KSG : :triangular +end + +_patch_ray_directions(lattice::LatticeType, patch::_LatticePatch) = + _lattice_directions(lattice)[patch.rays] """ - search_unweighted_gadgets(target_graph, target_boundary, loader; kwargs...) + check_crossing_frame(lattice, coordinates, pins, pin_rays) -Search for unweighted gadget replacements of `target_graph` by iterating over a `GraphLoader`. +Check the four geometric crossing-frame conditions. `pin_rays[i]` is the +outward lattice direction attached to `pins[i]`. The returned named tuple reports +G1 (strict hull interfaces), G2 (alternating channels), G3 (outward rays), and +G4 (clear pairwise non-adjacent exterior corridors). """ -function search_unweighted_gadgets( - target_graph::SimpleGraph{Int}, - target_boundary::Vector{Int}, - loader::GraphLoader; - prefilter::Bool=true, - limit::Union{Int,Nothing}=nothing, - max_results::Union{Int,Nothing}=nothing, +function check_crossing_frame( + lattice::LatticeType, + coordinates::Vector{_LatticeCoordinate}, + pins::Vector{_LatticeCoordinate}, + pin_rays::Vector{_LatticeCoordinate}, +) + length(pins) == 4 || throw(ArgumentError("a crossing frame requires four ordered pins")) + length(pin_rays) == 4 || throw(ArgumentError("a crossing frame requires four pin rays")) + length(unique(coordinates)) == length(coordinates) || + throw(ArgumentError("crossing-frame coordinates must be unique")) + length(unique(pins)) == 4 || + throw(ArgumentError("crossing-frame pins must be distinct")) + all(in(Set(coordinates)), pins) || + throw(ArgumentError("every crossing-frame pin must be present in coordinates")) + directions = _lattice_directions(lattice) + ray_indices = [_lattice_direction_index(directions, ray) for ray in pin_rays] + patch = _normalize_lattice_patch(lattice, coordinates, pins, ray_indices) + checks = _check_crossing_frame(lattice, patch) + return (G1=checks[1], G2=checks[2], G3=checks[3], G4=checks[4]) +end + +function _lattice_direction_index(directions, ray) + index = findfirst(==(ray), directions) + index === nothing && throw(ArgumentError("$ray is not a lattice direction")) + return index +end + +function _check_crossing_frame(lattice::LatticeType, patch::_LatticePatch) + length(patch.pins) == 4 || return (true, true, true, true) + directions = _patch_ray_directions(lattice, patch) + interfaces = [ + _lattice_step(lattice, pin, direction, 1) + for (pin, direction) in zip(patch.pins, directions) + ] + occupied_geometry = _geometry_coordinate.(Ref(lattice), patch.coordinates) + interface_geometry = _geometry_coordinate.(Ref(lattice), interfaces) + hull = _strict_convex_hull([occupied_geometry; interface_geometry]) + + g1 = length(unique(interface_geometry)) == 4 && all(in(hull), interface_geometry) + g2 = g1 && _interfaces_alternate(hull, interface_geometry) + g3 = _rays_point_outward(lattice, interfaces, directions) + g4 = _corridors_are_clear(lattice, patch.coordinates, patch.pins, interfaces, directions) + return (g1, g2, g3, g4) +end + +_canonical_coordinate(::Square, point::_LatticeCoordinate) = point +_canonical_coordinate(::Triangular, point::_LatticeCoordinate) = _offset_to_axial(point) +_from_canonical(::Square, point::_LatticeCoordinate) = point +_from_canonical(::Triangular, point::_LatticeCoordinate) = _axial_to_offset(point) +_lattice_distance(::Square, first, second) = max(abs(first[1] - second[1]), abs(first[2] - second[2])) +_lattice_distance(::Triangular, first, second) = max(abs(first[1] - second[1]), abs(first[2] - second[2]), abs(sum(first) - sum(second))) +_geometry_coordinate(::Square, point::_LatticeCoordinate) = point +function _geometry_coordinate(::Triangular, point::_LatticeCoordinate) + q, r = _offset_to_axial(point) + return (2q + r, r) +end + +_orientation(a, b, c) = + (b[1] - a[1]) * (c[2] - a[2]) - (b[2] - a[2]) * (c[1] - a[1]) + +function _strict_convex_hull(points::Vector{_LatticeCoordinate}) + sorted_points = sort!(unique(points)) + length(sorted_points) <= 2 && return sorted_points + lower = _LatticeCoordinate[] + for point in sorted_points + while length(lower) >= 2 && _orientation(lower[end-1], lower[end], point) <= 0 + pop!(lower) + end + push!(lower, point) + end + upper = _LatticeCoordinate[] + for point in Iterators.reverse(sorted_points) + while length(upper) >= 2 && _orientation(upper[end-1], upper[end], point) <= 0 + pop!(upper) + end + push!(upper, point) + end + return [lower[1:end-1]; upper[1:end-1]] +end + +function _interfaces_alternate( + hull::Vector{_LatticeCoordinate}, + interfaces::Vector{_LatticeCoordinate}, +) + labels = Dict(point => label for (label, point) in enumerate(interfaces)) + order = [labels[point] for point in hull if haskey(labels, point)] + length(order) == 4 || return false + return all(isodd(order[index]) != isodd(order[mod1(index + 1, 4)]) for index in 1:4) +end + +function _rays_point_outward( + lattice::LatticeType, + interfaces::Vector{_LatticeCoordinate}, + directions::Vector{_LatticeCoordinate}, +) + canonical_interfaces = _canonical_coordinate.(Ref(lattice), interfaces) + sum_q = sum(first, canonical_interfaces) + sum_r = sum(last, canonical_interfaces) + for (interface, direction) in zip(canonical_interfaces, directions) + out_q = 4interface[1] - sum_q + out_r = 4interface[2] - sum_r + if lattice isa Square + out_q * direction[1] + out_r * direction[2] > 0 || return false + else + out_x = 2out_q + out_r + direction_x = 2direction[1] + direction[2] + out_x * direction_x + 3out_r * direction[2] > 0 || return false + end + end + return true +end + +function _corridors_are_clear( + lattice::LatticeType, + coordinates::Vector{_LatticeCoordinate}, + pins::Vector{_LatticeCoordinate}, + interfaces::Vector{_LatticeCoordinate}, + directions::Vector{_LatticeCoordinate}, ) - total = isnothing(limit) ? length(loader) : min(length(loader), limit) - filter_fn = _make_unweighted_filter(target_graph, target_boundary; prefilter) - results = UnweightedGadget[] - @showprogress for key in Iterators.take(keys(loader), total) - result = filter_fn(loader[key], loader.layout[key], loader.pinset) - result === nothing && continue - push!(results, result) - max_results !== nothing && length(results) >= max_results && break - end - return results + occupied = _canonical_coordinate.(Ref(lattice), coordinates) + canonical_pins = _canonical_coordinate.(Ref(lattice), pins) + starts = _canonical_coordinate.(Ref(lattice), interfaces) + adjacency_offsets = [_LatticeCoordinate[(0, 0)]; _lattice_directions(lattice)] + + for index in eachindex(starts), site in occupied + site == canonical_pins[index] && continue + for offset in adjacency_offsets + _point_on_ray((site[1] + offset[1], site[2] + offset[2]), starts[index], directions[index]) && + return false + end + end + for first in 1:3, second in first+1:4, offset in adjacency_offsets + _rays_touch( + starts[first], directions[first], starts[second], directions[second], offset, + ) && return false + end + return true +end + +function _point_on_ray(point, start, direction) + displacement = (point[1] - start[1], point[2] - start[2]) + multiple = _direction_multiple(displacement, direction) + return multiple !== nothing && multiple >= 0 +end + +function _direction_multiple(displacement, direction) + if direction[1] != 0 + rem(displacement[1], direction[1]) == 0 || return nothing + multiple = div(displacement[1], direction[1]) + else + direction[2] != 0 || error("zero ray direction") + rem(displacement[2], direction[2]) == 0 || return nothing + multiple = div(displacement[2], direction[2]) + end + displacement == (multiple * direction[1], multiple * direction[2]) || return nothing + return multiple +end + +function _rays_touch(start1, direction1, start2, direction2, offset) + right = (start2[1] + offset[1] - start1[1], start2[2] + offset[2] - start1[2]) + determinant = direction1[2] * direction2[1] - direction1[1] * direction2[2] + if determinant != 0 + first_numerator = right[2] * direction2[1] - right[1] * direction2[2] + second_numerator = direction1[1] * right[2] - direction1[2] * right[1] + rem(first_numerator, determinant) == 0 || return false + rem(second_numerator, determinant) == 0 || return false + return div(first_numerator, determinant) >= 0 && div(second_numerator, determinant) >= 0 + end + + if direction1 == direction2 + return _direction_multiple(right, direction1) !== nothing + end + multiple = _direction_multiple(right, direction1) + return multiple !== nothing && multiple >= 0 end # ============================================================================ diff --git a/test/core/search.jl b/test/core/search.jl index c691f8e..a218bf7 100644 --- a/test/core/search.jl +++ b/test/core/search.jl @@ -313,7 +313,6 @@ end # Check that save file is created when results are found all_results = vcat(results...) if !isempty(all_results) && isfile(temp_file) - # Windows: avoid holding an open file handle during cleanup saved_data = JSON3.read(read(temp_file, String)) @test saved_data isa Union{Vector, JSON3.Array} @test length(saved_data) > 0 @@ -321,7 +320,7 @@ end finally # Clean up temp file - isfile(temp_file) && rm(temp_file; force=true, allow_delayed_delete=true) + isfile(temp_file) && rm(temp_file; force=true) end end diff --git a/test/core/unweighted_search.jl b/test/core/unweighted_search.jl index 3a7fa70..6527d44 100644 --- a/test/core/unweighted_search.jl +++ b/test/core/unweighted_search.jl @@ -1,128 +1,648 @@ using GadgetSearch +using Combinatorics using Graphs using Test -function _cross_graph() - g = SimpleGraph(4) - add_edge!(g, 1, 3) - add_edge!(g, 2, 4) - return g +function cross_graph() + graph = SimpleGraph(4) + add_edge!(graph, 1, 3) + add_edge!(graph, 2, 4) + return graph end -function _batoidea_graph() - g = SimpleGraph(11) - add_edge!(g, 1, 5); add_edge!(g, 1, 9) - add_edge!(g, 2, 5); add_edge!(g, 2, 6); add_edge!(g, 2, 7) - add_edge!(g, 3, 8) - add_edge!(g, 4, 9); add_edge!(g, 4, 10); add_edge!(g, 4, 11) - add_edge!(g, 5, 6); add_edge!(g, 5, 9); add_edge!(g, 5, 10) - add_edge!(g, 6, 7); add_edge!(g, 6, 9); add_edge!(g, 6, 10); add_edge!(g, 6, 11) - add_edge!(g, 7, 8); add_edge!(g, 7, 10); add_edge!(g, 7, 11) - add_edge!(g, 8, 11) - add_edge!(g, 9, 10) - add_edge!(g, 10, 11) - return g +function cross_edge_graph() + graph = SimpleGraph(6) + for edge in ((1,2), (1,5), (2,6), (3,5), (4,6)) + add_edge!(graph, edge...) + end + return graph end +@testset "Unweighted search" begin + @testset "selected connectivity matches enumeration" begin + lattice_host, _, _ = GadgetSearch._materialize_lattice_patch( + Triangular(), GadgetSearch._LatticePatch( + [(0,0),(1,0),(0,1),(1,1),(2,1)], [(0,0)], [1], + ), + ) + for host in (path_graph(5), cycle_graph(5), lattice_host) + adjacent = [collect(neighbors(host, vertex)) for vertex in vertices(host)] + for mask in 0:(1 << nv(host))-1 + !iszero(mask & 1) || continue + chosen = [vertex for vertex in vertices(host) + if !iszero(mask & (1 << (vertex - 1)))] + cnf = GadgetSearch._SatCnf() + selected = [GadgetSearch._sat_variable!(cnf) for _ in vertices(host)] + GadgetSearch._add_selected_connectivity!( + cnf, selected, adjacent, 1, length(chosen), + ) + for vertex in vertices(host) + GadgetSearch._sat_clause!( + cnf, vertex in chosen ? selected[vertex] : -selected[vertex], + ) + end + actual = GadgetSearch._next_selected_assignment!( + GadgetSearch._new_sat_solver(cnf), selected, + ) !== nothing + @test actual == is_connected(induced_subgraph(host, chosen)[1]) + end + end + end -function _edge_graph() - g = SimpleGraph(2) - add_edge!(g, 1, 2) - return g -end + @testset "fixed verifier" begin + reduced = calculate_reduced_alpha_tensor(cross_graph(), [1, 2, 3, 4]) + @test GadgetSearch.inf_mask(reduced) == BigInt(60576) + @test is_diff_by_constant(reduced .+ 3, reduced) == (true, 3.0) + end -function _connected_graph() - g = SimpleGraph(4) - add_edge!(g, 1, 3) - add_edge!(g, 1, 4) - add_edge!(g, 2, 4) - add_edge!(g, 3, 4) - return g -end + @testset "crossing frame" begin + square = [(0, 0), (-1, 0), (0, 1), (1, 0), (0, -1)] + pins = square[2:5] + rays = [(-1, 0), (0, 1), (1, 0), (0, -1)] + @test all(check_crossing_frame(Square(), square, pins, rays)) + @test !check_crossing_frame(Square(), [square; (-2, 1)], pins, rays).G4 + @test_throws ArgumentError check_crossing_frame( + Square(), [(0, 0)], pins, rays, + ) + @test_throws ArgumentError check_crossing_frame( + Square(), [square; square[1]], pins, rays, + ) + @test_throws ArgumentError check_crossing_frame( + Square(), square, [pins[1], pins[1], pins[3], pins[4]], rays, + ) + triangular = [(0, 0), (-1, 0), (0, 1), (1, 0), (-1, -1)] + triangular_pins = triangular[2:5] + @test all(check_crossing_frame( + Triangular(), triangular, triangular_pins, rays, + )) + end -function _isolated_graph() - g = SimpleGraph(3) - add_edge!(g, 1, 2) - return g -end + @testset "pruned frame generation preserves geometric candidates" begin + lattice = Triangular() + pins = [(2,0),(2,3),(0,1),(1,0)] + ray_indices = eachindex(GadgetSearch._lattice_directions(lattice)) + brute_force = Set{Tuple}() + for order in permutations(1:4) + ordered_pins = pins[collect(order)] + for rays in Iterators.product(ntuple(_ -> ray_indices, 4)...) + all(GadgetSearch._check_crossing_frame( + lattice, + GadgetSearch._LatticePatch( + ordered_pins, ordered_pins, collect(rays), + ), + )) || continue + push!(brute_force, (Tuple(ordered_pins), rays)) + end + end -function _to_g6(g) - return graph_to_g6(g) -end + pruned = Set{Tuple}() + ray_options = [ + [ray for ray in ray_indices if GadgetSearch._single_pin_corridor_clear( + lattice, pins, index, ray, + )] for index in eachindex(pins) + ] + for rays_tuple in Iterators.product(ray_options...) + GadgetSearch._pin_rays_are_compatible( + lattice, pins, rays_tuple, + ) || continue + physical_rays = collect(rays_tuple) + checks = GadgetSearch._check_crossing_frame( + lattice, GadgetSearch._LatticePatch(pins, pins, physical_rays), + ) + checks[1] && checks[3] && checks[4] || continue + for order_tuple in permutations(1:4) + order = collect(order_tuple) + ordered_pins = pins[order] + ordered_rays = physical_rays[order] + GadgetSearch._labels_alternate( + lattice, ordered_pins, ordered_rays, + ) || continue + push!(pruned, (Tuple(ordered_pins), Tuple(ordered_rays))) + end + end + @test pruned == brute_force -@testset "Unweighted Search" begin - @testset "search_unweighted_gadgets: basic" begin - cross = _cross_graph() - batoidea = _batoidea_graph() - loader = GraphLoader( - GraphDataset([_to_g6(cross), _to_g6(batoidea)]), - pinset=[1, 2, 3, 4], - ) - results = search_unweighted_gadgets(cross, [1, 2, 3, 4], loader) - @test results isa Vector{UnweightedGadget} - @test any(r -> r.constant_offset == 0.0, results) - @test any(r -> r.constant_offset == 2.0, results) - @test all(r -> r.pattern_graph == cross, results) - @test !hasproperty(UnweightedGadget, :target_index) - end + clique_frames = Set{Tuple}() + expected_prefixed_frames = Set{Tuple}() + clique_evaluated = Ref(0) + candidates = GadgetSearch._crossing_port_candidates( + lattice, [(column, row) for column in 0:2 for row in 0:2], + ) + GadgetSearch._foreach_crossing_frame_clique( + lattice, (3,3), clique_evaluated, typemax(Int); min_allowed=4, + ) do frame, _, _ + key = GadgetSearch._canonical_crossing_frame_key(lattice, frame) + push!(clique_frames, key) + indices = sort([findfirst(==(port), candidates) + for port in zip(frame.pins, frame.rays)]) + indices[1:2] == [1, 9] && push!(expected_prefixed_frames, key) + return nothing + end + sharded_frames = Set{Tuple}() + for shard_index in 0:2 + evaluated = Ref(0) + GadgetSearch._foreach_crossing_frame_clique( + lattice, (3,3), evaluated, typemax(Int); + min_allowed=4, shard_index, shard_count=3, + ) do frame, _, _ + push!(sharded_frames, + GadgetSearch._canonical_crossing_frame_key(lattice, frame)) + return nothing + end + end + @test sharded_frames == clique_frames - @testset "search_unweighted_gadgets: limit and max_results" begin - cross = _cross_graph() - batoidea = _batoidea_graph() - loader = GraphLoader( - GraphDataset([_to_g6(cross), _to_g6(batoidea)]), - pinset=[1, 2, 3, 4], - ) - limited = search_unweighted_gadgets(cross, [1, 2, 3, 4], loader; limit=1) - @test length(limited) == 1 - @test limited[1].constant_offset == 0.0 - capped = search_unweighted_gadgets(cross, [1, 2, 3, 4], loader; max_results=1) - @test length(capped) == 1 - end + prefixed_frames = Set{Tuple}() + prefixed_evaluated = Ref(0) + GadgetSearch._foreach_crossing_frame_clique( + lattice, (3,3), prefixed_evaluated, typemax(Int); + min_allowed=4, port_prefix=(1, 9), + ) do frame, _, _ + push!(prefixed_frames, + GadgetSearch._canonical_crossing_frame_key(lattice, frame)) + return nothing + end + @test prefixed_frames == expected_prefixed_frames + cross23_window = [(column, row) + for column in 0:7 for row in 0:7] + @test GadgetSearch._crossing_port_prefix( + Triangular(), cross23_window, 5117, + ) == (21, 177) + + for (candidate_pins, candidate_rays, side) in ( + (pins, [6,2,4,5], 4), + ([(5,1),(0,3),(3,5),(6,6)], [6,3,3,2], 8), + ) + window = [(column, row) + for column in 0:side-1 for row in 0:side-1] + expected = sort([site for site in window + if site in candidate_pins || begin + checks = GadgetSearch._check_crossing_frame( + lattice, GadgetSearch._LatticePatch( + [candidate_pins; site], candidate_pins, candidate_rays, + ), + ) + checks[1] && checks[3] && checks[4] + end]) + @test GadgetSearch._allowed_crossing_sites( + lattice, window, candidate_pins, candidate_rays, 0, + ) == expected + end - @testset "search_unweighted_gadgets: prefilter rejects disconnected pin coverage" begin - loader = GraphLoader(GraphDataset([_to_g6(_cross_graph())]), pinset=[1, 3]) - edge = _edge_graph() - results_on = search_unweighted_gadgets(edge, [1, 2], loader; prefilter=true) - results_off = search_unweighted_gadgets(edge, [1, 2], loader; prefilter=false) - @test isempty(results_on) - @test length(results_off) == 1 end - @testset "UnweightedGadget has no target_index" begin - @test !(:target_index in fieldnames(UnweightedGadget)) + @testset "fixed SAT positive control" begin + coordinates = [ + (0,3),(6,6),(1,3),(2,1),(2,2),(2,3),(3,1),(3,3), + (3,4),(3,5),(4,1),(4,2),(4,4),(4,5),(5,1),(5,2), + (5,3),(5,4),(5,5),(6,2),(6,3),(6,5),(7,4), + ] + frame = ( + pins=[(5,1),(0,3),(3,5),(6,6)], + rays=[6,3,3,2], + allowed=coordinates, + ) + target_reduced = vec(calculate_reduced_alpha_tensor( + cross_graph(), [1,2,3,4], + )) + analysis = GadgetSearch._solve_fixed_crossing_sat( + target_reduced, Triangular(), frame, 23, 7, + ) + @test analysis !== nothing + @test analysis.offset == 7 + context = GadgetSearch._prepare_sat_frame(Triangular(), frame) + filter_solver, filter_selected = GadgetSearch._lower_state_filter_problem( + target_reduced, context, 23, 7, + first(GadgetSearch._essential_lower_states(target_reduced)), + ) + @test GadgetSearch._next_selected_assignment!( + filter_solver, filter_selected, + ) !== nothing + @test nv(analysis.graph) == 23 + @test is_gadget_replacement( + cross_graph(), analysis.graph, [1,2,3,4], analysis.boundary, + ) == (true, 7.0) + @test all(check_crossing_frame( + Triangular(), analysis.patch.coordinates, analysis.patch.pins, + GadgetSearch._patch_ray_directions(Triangular(), analysis.patch), + )) + + extended_sites = [coordinates; (6,0); (3,6)] + extended_pins = [(6,0), (0,3), (3,6), (6,6)] + extended_analysis = GadgetSearch._analyze_crossing_candidate( + target_reduced, Triangular(), extended_sites, extended_pins, + frame.rays, + ) + @test extended_analysis.solved + @test extended_analysis.offset == 8 + extended_gadget = GadgetSearch._unweighted_gadget( + cross_graph(), Triangular(), extended_analysis, + ) + optimized = optimize_unweighted_gadget( + extended_gadget, [1,2,3,4]; + min_vertices=23, max_sat_evaluations=0, + ) + @test nv(optimized.gadget.replacement_graph) == 23 + @test optimized.gadget.constant_offset == 7 + @test only(optimized.steps).rule == :opposite_leaf_pin_contraction + @test nv(only(optimized.steps).before.replacement_graph) == 25 + @test nv(only(optimized.steps).after.replacement_graph) == 23 + @test optimized.termination_reason == :minimum_vertices + @test is_gadget_replacement( + cross_graph(), optimized.gadget.replacement_graph, [1,2,3,4], + optimized.gadget.boundary_vertices, + ) == (true, 7.0) + + directions = GadgetSearch._lattice_directions(Triangular()) + tail_middle = GadgetSearch._lattice_step( + Triangular(), frame.pins[1], directions[frame.rays[1]], 1, + ) + tail_pin = GadgetSearch._lattice_step( + Triangular(), frame.pins[1], directions[frame.rays[1]], 2, + ) + tail_analysis = GadgetSearch._analyze_crossing_candidate( + target_reduced, Triangular(), + [coordinates; tail_middle; tail_pin], + [tail_pin; frame.pins[2:4]], frame.rays, + ) + @test tail_analysis.solved + tail_gadget = GadgetSearch._unweighted_gadget( + cross_graph(), Triangular(), tail_analysis, + ) + tail_optimized = optimize_unweighted_gadget( + tail_gadget, [1,2,3,4]; + min_vertices=23, max_sat_evaluations=0, + ) + @test only(tail_optimized.steps).rule == + :even_boundary_tail_contraction + @test nv(tail_optimized.gadget.replacement_graph) == 23 + + rewrite_start_coordinates = [ + (4,7),(2,0),(6,0),(7,3),(1,3),(2,4),(3,6),(3,7), + (2,2),(3,4),(3,5),(2,1),(3,3),(4,4),(4,5),(3,0), + (4,2),(5,4),(5,5),(4,0),(4,1),(5,2),(5,3),(6,5), + (5,1),(6,2),(6,3),(7,4), + ] + rewrite_start_analysis = GadgetSearch._analyze_crossing_candidate( + target_reduced, Triangular(), rewrite_start_coordinates, + [(4,7),(2,0),(6,0),(7,3)], [2,5,5,1], + ) + @test rewrite_start_analysis.solved + @test rewrite_start_analysis.offset == 10 + rewrite_start = GadgetSearch._unweighted_gadget( + cross_graph(), Triangular(), rewrite_start_analysis, + ) + rewritten = optimize_unweighted_gadget( + rewrite_start, [1,2,3,4]; + min_vertices=23, max_sat_evaluations=8, + max_sat_conflicts=100_000, host_radius=1, + ) + @test rewritten.termination_reason == :minimum_vertices + @test rewritten.sat_evaluations <= 8 + @test rewritten.unresolved_sat_evaluations < rewritten.sat_evaluations + @test [step.rule for step in rewritten.steps] == [ + :even_boundary_tail_contraction, + :frame_rewrite_resynthesis, + :opposite_leaf_pin_contraction, + ] + @test [ + nv(step.before.replacement_graph) => nv(step.after.replacement_graph) + for step in rewritten.steps + ] == [28 => 26, 26 => 25, 25 => 23] + @test is_gadget_replacement( + cross_graph(), rewritten.gadget.replacement_graph, [1,2,3,4], + rewritten.gadget.boundary_vertices, + ) == (true, 7.0) + + rewrite_budget = optimize_unweighted_gadget( + GadgetSearch._unweighted_gadget( + cross_graph(), Triangular(), analysis, + ), + [1,2,3,4]; + min_vertices=22, max_sat_evaluations=1, host_radius=0, + ) + @test isempty(rewrite_budget.steps) + @test rewrite_budget.sat_evaluations == 1 + @test rewrite_budget.termination_reason == :sat_budget + + rewrite_unknown = optimize_unweighted_gadget( + GadgetSearch._unweighted_gadget( + cross_graph(), Triangular(), analysis, + ), + [1,2,3,4]; + min_vertices=22, max_sat_evaluations=10_000, + max_sat_conflicts=0, host_radius=0, + ) + @test isempty(rewrite_unknown.steps) + @test rewrite_unknown.sat_evaluations < 10_000 + @test rewrite_unknown.unresolved_sat_evaluations > 0 + @test rewrite_unknown.termination_reason == :sat_unknown + canonical_coordinates = [ + (0,-5),(-1,-4),(-1,-3),(1,-5),(0,-3),(3,-7),(2,-6), + (2,-5),(1,-4),(1,-3),(0,-2),(0,-1),(3,-6),(2,-4), + (1,-2),(0,0),(5,-7),(4,-6),(3,-4),(3,-3),(2,-2), + (6,-7),(4,-4),(4,-3),(6,-6),(6,-5),(5,-4),(4,-2), + ] + joint_frame = [ + ((0,0),(1,0)), ((-1,-3),(-1,1)), + ((3,-7),(0,-1)), ((4,-2),(1,0)), + ] + joint_solution = GadgetSearch._from_canonical.( + Ref(Triangular()), canonical_coordinates, + ) + joint_pins = GadgetSearch._from_canonical.( + Ref(Triangular()), first.(joint_frame), + ) + directions = GadgetSearch._lattice_directions(Triangular()) + joint_rays = [findfirst(==(ray), directions) for ray in last.(joint_frame)] + joint_analysis = GadgetSearch._analyze_crossing_candidate( + target_reduced, Triangular(), joint_solution, joint_pins, joint_rays, + ) + @test joint_analysis.solved + @test joint_analysis.offset == 10 + joint_cnf, joint_selected, joint_choices, joint_coordinates = + GadgetSearch._joint_crossing_sat_cnf( + target_reduced, Triangular(), (8,8), 28, 10, + ) + selected_coordinates = Set(joint_solution) + for (vertex, coordinate) in enumerate(joint_coordinates) + GadgetSearch._sat_clause!( + joint_cnf, coordinate in selected_coordinates ? + joint_selected[vertex] : -joint_selected[vertex], + ) + end + for label in 1:4, + (vertex, direction, variable) in joint_choices[label] + chosen = joint_coordinates[vertex] == joint_pins[label] && + direction == joint_rays[label] + GadgetSearch._sat_clause!( + joint_cnf, chosen ? variable : -variable, + ) + end + status, _ = GadgetSearch._next_selected_assignment_limited!( + GadgetSearch._new_sat_solver(joint_cnf), joint_selected, 100_000, + ) + @test status == :sat + + cross23_sites = [ + (0,0),(3,-7),(0,-1),(-2,-1),(-1,-1),(0,-2),(-2,-2), + (0,-3),(1,-3),(2,-4),(-2,-3),(-1,-3),(1,-4),(2,-5), + (-2,-4),(-1,-4),(0,-5),(1,-5),(2,-6),(-1,-5),(0,-6), + (2,-7),(1,-7), + ] + cross23_pins = [(0,0),(-2,-4),(3,-7),(2,-4)] + cross23_rays = [1,4,6,1] + cross23_solution = GadgetSearch._from_canonical.( + Ref(Triangular()), cross23_sites, + ) + cross23_boundary = GadgetSearch._from_canonical.( + Ref(Triangular()), cross23_pins, + ) + cross23_cnf, cross23_selected, cross23_choices, cross23_coordinates = + GadgetSearch._joint_crossing_sat_cnf( + target_reduced, Triangular(), (8,8), 23, 7; + canonical_shift=(-1,0), + ) + selected_coordinates = Set(cross23_solution) + for (vertex, coordinate) in enumerate(cross23_coordinates) + GadgetSearch._sat_clause!( + cross23_cnf, coordinate in selected_coordinates ? + cross23_selected[vertex] : -cross23_selected[vertex], + ) + end + for label in 1:4, + (vertex, direction, variable) in cross23_choices[label] + chosen = cross23_coordinates[vertex] == cross23_boundary[label] && + direction == cross23_rays[label] + GadgetSearch._sat_clause!( + cross23_cnf, chosen ? variable : -variable, + ) + end + status, _ = GadgetSearch._next_selected_assignment_limited!( + GadgetSearch._new_sat_solver(cross23_cnf), cross23_selected, + 100_000, + ) + @test status == :sat + cross23_analysis = GadgetSearch._analyze_crossing_candidate( + target_reduced, Triangular(), cross23_solution, + cross23_boundary, cross23_rays, + ) + @test cross23_analysis.solved + @test cross23_analysis.offset == 7 + shifted_target = map(value -> isfinite(value) ? value + 8 : value, target_reduced) + negative_offset = GadgetSearch._solve_fixed_crossing_sat( + shifted_target, Triangular(), frame, 23, -1, + ) + @test negative_offset !== nothing + @test negative_offset.offset == -1 + + ksg_coordinates = [(0, 0), (-1, 0), (0, 1), (1, 0), (0, -1)] + ksg_frame = ( + pins=ksg_coordinates[2:5], + rays=[4, 7, 5, 2], + allowed=ksg_coordinates, + ) + ksg_patch = GadgetSearch._LatticePatch( + ksg_coordinates, ksg_frame.pins, ksg_frame.rays, + ) + ksg_graph, ksg_boundary, _ = GadgetSearch._materialize_lattice_patch( + Square(), ksg_patch, + ) + ksg_target = vec(calculate_reduced_alpha_tensor(ksg_graph, ksg_boundary)) + ksg_analysis = GadgetSearch._solve_fixed_crossing_sat( + ksg_target, Square(), ksg_frame, 5, 0, + ) + @test ksg_analysis !== nothing + @test nv(ksg_analysis.graph) == 5 + @test all(check_crossing_frame( + Square(), ksg_analysis.patch.coordinates, ksg_analysis.patch.pins, + GadgetSearch._patch_ray_directions(Square(), ksg_analysis.patch), + )) + @test isnothing(GadgetSearch._solve_fixed_crossing_sat( + ksg_target, Square(), ksg_frame, 6, 0, + )) end - @testset "inf_mask (internal)" begin - @test GadgetSearch.inf_mask([0.0, -Inf, 3.0, -Inf]) == BigInt(10) - @test GadgetSearch.inf_mask(fill(-Inf, 4)) == BigInt(15) - reduced = calculate_reduced_alpha_tensor(_cross_graph(), [1, 2, 3, 4]) - @test GadgetSearch.inf_mask(reduced) == BigInt(60576) + @testset "joint frame and occupancy SAT" begin + target = complete_graph(4) + target_reduced = vec(calculate_reduced_alpha_tensor( + target, collect(1:4), + )) + public_joint = search_unweighted_gadget_joint( + target, collect(1:4), Square(); + window_shape=(2,2), atom_count=4, offset=0, seconds=30, + ) + @test public_joint !== nothing + @test nv(public_joint.replacement_graph) == 4 + @test first(is_gadget_replacement( + target, public_joint.replacement_graph, collect(1:4), + public_joint.boundary_vertices, + )) + @test all(check_crossing_frame( + Square(), public_joint.lattice_coordinates, + public_joint.lattice_coordinates[public_joint.boundary_vertices], + public_joint.pin_rays, + )) + cross_edge_target = cross_edge_graph() + cross_edge = search_unweighted_gadget_joint( + cross_edge_target, collect(1:4), Triangular(); + window_shape=(3,4), atom_count=9, offset=1, seconds=30, + canonical_shift=(-1,1), first_ray=1, + ) + @test cross_edge !== nothing + @test first(is_gadget_replacement( + cross_edge_target, cross_edge.replacement_graph, collect(1:4), + cross_edge.boundary_vertices, + )) + optimized_cross_edge = optimize_unweighted_gadget( + cross_edge, collect(1:4); + min_vertices=4, max_sat_evaluations=256, + max_sat_conflicts=100_000, host_radius=1, + ) + @test nv(optimized_cross_edge.gadget.replacement_graph) == 9 + @test isempty(optimized_cross_edge.steps) + @test optimized_cross_edge.sat_evaluations == 256 + @test optimized_cross_edge.termination_reason == :sat_budget + rotated_cross_edge = search_unweighted_gadget_joint( + cross_edge_target, collect(1:4), Triangular(); + window_shape=(3,4), atom_count=9, offset=1, seconds=30, + canonical_shift=(0,0), first_ray=2, + ) + @test rotated_cross_edge !== nothing + @test first(is_gadget_replacement( + cross_edge_target, rotated_cross_edge.replacement_graph, + collect(1:4), rotated_cross_edge.boundary_vertices, + )) + @test isnothing(GadgetSearch._solve_joint_crossing_sat( + target_reduced, Square(), (2,2), 4, 1, + )) + + crossing_points = [(0,0), (0,2), (2,0), (2,2)] + for chosen in Iterators.product(ntuple(_ -> 1:4, 4)...) + points = crossing_points[collect(chosen)] + hull = GadgetSearch._strict_convex_hull(points) + expected = length(unique(points)) == 4 && length(hull) == 4 && + GadgetSearch._interfaces_alternate(hull, points) + actual = GadgetSearch._interfaces_form_alternating_quadrilateral( + points, + ) + @test actual == expected + end end - @testset "pins_prefilter (internal)" begin - connected = _connected_graph() - disconnected = _cross_graph() - isolated = _isolated_graph() - @test GadgetSearch.pins_prefilter(connected, [1]) - @test GadgetSearch.pins_prefilter(disconnected, [1, 2]) - @test !GadgetSearch.pins_prefilter(disconnected, [1]) - @test !GadgetSearch.pins_prefilter(isolated, [1]) - @test GadgetSearch.pins_prefilter(isolated, [1, 3]) - @test_throws ErrorException GadgetSearch.pins_prefilter(connected, [1, 1]) - @test_throws ErrorException GadgetSearch.pins_prefilter(connected, [0]) + @testset "SAT model enumeration" begin + cnf = GadgetSearch._SatCnf() + foreach(_ -> GadgetSearch._sat_variable!(cnf), 1:3) + solver = GadgetSearch._new_sat_solver(cnf) + assignments = Set{Tuple{Bool, Bool}}() + while true + assignment = GadgetSearch._next_selected_assignment!(solver, [1, 2]) + assignment === nothing && break + push!(assignments, Tuple(assignment)) + end + @test assignments == Set([ + (false, false), (false, true), (true, false), (true, true), + ]) + sat_cnf = GadgetSearch._SatCnf() + variable = GadgetSearch._sat_variable!(sat_cnf) + GadgetSearch._sat_clause!(sat_cnf, variable) + status, assignment = GadgetSearch._next_selected_assignment_limited!( + GadgetSearch._new_sat_solver(sat_cnf), [variable], 10, + ) + @test status == :sat + @test assignment == [true] + + unsat_cnf = GadgetSearch._SatCnf() + variable = GadgetSearch._sat_variable!(unsat_cnf) + GadgetSearch._sat_clause!(unsat_cnf, variable) + GadgetSearch._sat_clause!(unsat_cnf, -variable) + status, assignment = GadgetSearch._next_selected_assignment_limited!( + GadgetSearch._new_sat_solver(unsat_cnf), [variable], 10, + ) + @test status == :unsat + @test assignment === nothing end - @testset "Triangular UDG Integration" begin - path = tempname() * ".g6" - try - generate_full_grid_udg(Triangular(), 1, 1; path=path) - loader = GraphLoader(path; pinset=[1, 2, 3, 4]) - target = loader[1] - results = search_unweighted_gadgets(target, [1, 2, 3, 4], loader; limit=1, max_results=1) - @test length(results) == 1 - @test results[1].constant_offset == 0.0 - finally - isfile(path) && rm(path) + @testset "public bounded search" begin + @test GadgetSearch._crossing_window_shapes(8, 23) == [ + (8,8), (8,7), (7,8), (8,6), (6,8), + (8,5), (5,8), (8,4), (4,8), (8,3), (3,8), + ] + @test_throws ArgumentError search_unweighted_gadgets( + path_graph(3), [1,2,3], Triangular(), + ) + @test_throws ArgumentError search_unweighted_gadgets( + path_graph(4), [1,2,3,4]; window_side=1, + ) + @test_throws ArgumentError search_unweighted_gadgets( + path_graph(4), [1,2,3,4]; checkpoint_interval=0, + ) + ksg = search_unweighted_gadgets( + cross_graph(), [1,2,3,4], Square(); + min_vertices=4, max_vertices=4, max_evaluations=1, + ) + @test ksg.lattice == :KSG + @test ksg.evaluated <= 1 + @test ksg.termination_reason in (:budget, :frame_budget) + + exhausted = search_unweighted_gadgets( + cross_graph(), [1,2,3,4], Triangular(); + min_vertices=4, max_vertices=4, max_evaluations=1, + ) + @test exhausted.evaluated == 1 + @test isempty(exhausted.gadgets) + @test exhausted.termination_reason in (:budget, :frame_budget) + + ksg_target = complete_graph(4) + solved = search_unweighted_gadgets( + ksg_target, collect(1:4), Square(); + min_vertices=4, max_vertices=4, max_evaluations=50, max_results=1, + max_frame_evaluations=100_000, window_side=2, + ) + @test solved.termination_reason == :solution + @test length(solved.gadgets) == 1 + @test is_gadget_replacement( + ksg_target, solved.gadgets[1].replacement_graph, + collect(1:4), solved.gadgets[1].boundary_vertices, + )[1] + + fully_enumerated = search_unweighted_gadgets( + ksg_target, collect(1:4), Square(); + min_vertices=4, max_vertices=4, max_evaluations=50, max_results=2, + max_frame_evaluations=100_000, window_side=2, + ) + @test fully_enumerated.termination_reason == :search_space_exhausted + @test length(fully_enumerated.gadgets) == 1 + @test fully_enumerated.evaluated > solved.evaluated + + mktempdir() do directory + checkpoint_path = joinpath(directory, "search.checkpoint") + interrupted = search_unweighted_gadgets( + ksg_target, collect(1:4), Square(); + min_vertices=4, max_vertices=4, max_evaluations=3, + max_results=2, max_frame_evaluations=100_000, + window_side=2, checkpoint_path, checkpoint_interval=1, + ) + @test interrupted.termination_reason == :budget + @test isfile(checkpoint_path) + checkpoint = read_unweighted_search_checkpoint(checkpoint_path) + @test checkpoint.window_shape == (2, 2) + @test checkpoint.frame_cursor >= 1 + @test checkpoint.order_cursor >= 1 + @test checkpoint.offset_cursor >= 0 + resumed = search_unweighted_gadgets( + ksg_target, collect(1:4), Square(); + min_vertices=4, max_vertices=4, max_evaluations=20, + max_results=2, max_frame_evaluations=100_000, + window_side=2, checkpoint_path, checkpoint_interval=1, + ) + @test resumed.termination_reason == :search_space_exhausted + @test length(resumed.gadgets) == 1 end + point = (3, 4) + @test GadgetSearch._from_canonical(Square(), point) == point + canonical = GadgetSearch._canonical_coordinate(Triangular(), point) + @test GadgetSearch._from_canonical(Triangular(), canonical) == point + @test GadgetSearch._lattice_distance(Square(), (0, 0), (2, -1)) == 2 + @test GadgetSearch._lattice_distance(Triangular(), (0, 0), (2, -1)) == 2 end end