diff --git a/CHANGELOG.md b/CHANGELOG.md index 8d2f60a..1875e6b 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,29 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] +### Added + +- `MinimizationComplexityWarning`, raised by `minimize` once the primes are + known but before the exponential cover search begins, so a run that will take + a long time says so rather than appearing to hang. The run still completes and + the answer is still exact; silence it with `warnings.simplefilter` or disable + it with `minimize(..., complexity_guard=False)`. The threshold is the prime + count, chosen from measurement rather than intuition. +- `benchmarks/profile_phases.py`, which times truth-table construction, prime + generation, chart construction and cover solving separately across cases, + conditions, sufficient share and remainder share. Remainder-heavy problems are + dominated by prime generation; dense on-sets by the cover search. +- `build_chart` accepts pre-computed `primes`, so callers that already generated + them do not pay for it twice. + +### Changed + +- Prime generation now works on integer bitmasks rather than tuples of optional + bits: **141× faster** at eight conditions (1.41 s to 0.010 s), 83× end to end. + The algorithm is unchanged and every exactness and R-parity test passes + unchanged. A ten-condition parsimonious solution now takes 0.45 s rather than + roughly 36 s. + ## [0.2.0] — 2026-08-11 Archived on Zenodo: [10.5281/zenodo.21887472](https://doi.org/10.5281/zenodo.21887472) diff --git a/benchmarks/profile_phases.py b/benchmarks/profile_phases.py new file mode 100644 index 0000000..2ccdd54 --- /dev/null +++ b/benchmarks/profile_phases.py @@ -0,0 +1,307 @@ +"""Phase-level benchmarks for the QCA pipeline. + +Measures each stage separately so that optimisation follows evidence rather +than intuition: + +1. truth-table construction +2. prime-implicant generation +3. chart construction +4. exact cover solving +5. end-to-end minimisation + +Four dimensions are varied independently, because they do not cost the same: +the number of cases, the number of conditions, the number of sufficient +configurations, and the number of logical remainders. + +Usage +----- + python benchmarks/profile_phases.py # the default sweep + python benchmarks/profile_phases.py --max-width 9 # push further + python benchmarks/profile_phases.py --markdown # table for the docs + +Runtime grows sharply with the number of conditions, so the defaults are +chosen to finish in seconds and are safe to run in CI. +""" + +from __future__ import annotations + +import argparse +from dataclasses import dataclass +from random import Random +from time import perf_counter + +import pandas as pd + +from setqca import build_truth_table +from setqca.minimize import build_chart, prime_implicants +from setqca.minimize.qmc import solve_minimum_cover + + +@dataclass(frozen=True, slots=True) +class Timing: + """One measured configuration.""" + + label: str + cases: int + conditions: int + positives: int + remainders: int + primes: int + truth_table_s: float + primes_s: float + chart_s: float + cover_s: float + + @property + def total_s(self) -> float: + """Return the summed time across the measured phases.""" + return self.truth_table_s + self.primes_s + self.chart_s + self.cover_s + + @property + def dominant(self) -> str: + """Return the phase taking the largest share.""" + phases = { + "truth table": self.truth_table_s, + "primes": self.primes_s, + "chart": self.chart_s, + "cover": self.cover_s, + } + return max(phases, key=lambda name: phases[name]) + + +def _fuzzy_frame(rng: Random, cases: int, conditions: int) -> pd.DataFrame: + """Build calibrated data whose outcome actually depends on the conditions. + + An outcome drawn independently of the conditions makes almost every + observed row consistent, so the minimiser collapses to the tautology and + the benchmark measures nothing. Here the outcome follows two overlapping + conjunctions plus noise, which is the shape QCA data normally has. + """ + columns = { + f"C{index}": [round(rng.uniform(0.02, 0.98), 3) for _ in range(cases)] + for index in range(conditions) + } + outcome: list[float] = [] + for case in range(cases): + memberships = [columns[f"C{index}"][case] for index in range(conditions)] + first = min(memberships[: max(2, conditions // 2)]) + second = min(memberships[-2:]) if conditions >= 2 else memberships[0] + noisy = max(first, second) + rng.uniform(-0.25, 0.25) + outcome.append(round(min(0.98, max(0.02, noisy)), 3)) + columns["Y"] = outcome + frame = pd.DataFrame(columns) + # Nudge anything that landed on 0.5, which the truth table rejects. + return frame.mask(frame == 0.5, 0.51) + + +def measure_truth_table(rng: Random, cases: int, conditions: int) -> tuple[float, object]: + """Time truth-table construction for a random calibrated dataset.""" + frame = _fuzzy_frame(rng, cases, conditions) + names = [f"C{index}" for index in range(conditions)] + start = perf_counter() + table = build_truth_table(frame, outcome="Y", conditions=names, inclusion_cutoff=0.75) + return perf_counter() - start, table + + +def measure_minimisation( + on_set: set[int], dont_cares: set[int], width: int +) -> tuple[float, float, float, int]: + """Time prime generation, chart construction and cover solving separately.""" + start = perf_counter() + primes = prime_implicants(on_set, dont_cares, width) + primes_s = perf_counter() - start + + start = perf_counter() + chart = build_chart(on_set, dont_cares=dont_cares, width=width) + chart_s = perf_counter() - start + + covered = [frozenset(m for m in on_set if prime.covers(m)) for prime in primes] + literals = [prime.literals for prime in primes] + start = perf_counter() + solve_minimum_cover(covered, literals, on_set, max_solutions=32) + cover_s = perf_counter() - start + + return primes_s, chart_s, cover_s, len(chart.primes) + + +def sweep_conditions(rng: Random, max_width: int, cases: int) -> list[Timing]: + """Vary the number of conditions, holding the case count fixed.""" + results: list[Timing] = [] + for width in range(3, max_width + 1): + truth_table_s, table = measure_truth_table(rng, cases, width) + on_set = table.positive_minterms # type: ignore[attr-defined] + remainders = table.remainder_minterms # type: ignore[attr-defined] + if not on_set: + continue + primes_s, chart_s, cover_s, primes = measure_minimisation(on_set, remainders, width) + results.append( + Timing( + label=f"conditions={width}", + cases=cases, + conditions=width, + positives=len(on_set), + remainders=len(remainders), + primes=primes, + truth_table_s=truth_table_s, + primes_s=primes_s, + chart_s=chart_s, + cover_s=cover_s, + ) + ) + return results + + +def sweep_cases(rng: Random, conditions: int, counts: tuple[int, ...]) -> list[Timing]: + """Vary the number of cases, holding the property space fixed.""" + results: list[Timing] = [] + for cases in counts: + truth_table_s, table = measure_truth_table(rng, cases, conditions) + on_set = table.positive_minterms # type: ignore[attr-defined] + remainders = table.remainder_minterms # type: ignore[attr-defined] + if not on_set: + continue + primes_s, chart_s, cover_s, primes = measure_minimisation(on_set, remainders, conditions) + results.append( + Timing( + label=f"cases={cases}", + cases=cases, + conditions=conditions, + positives=len(on_set), + remainders=len(remainders), + primes=primes, + truth_table_s=truth_table_s, + primes_s=primes_s, + chart_s=chart_s, + cover_s=cover_s, + ) + ) + return results + + +def sweep_density(rng: Random, width: int, densities: tuple[float, ...]) -> list[Timing]: + """Vary how many configurations are sufficient, at a fixed width.""" + results: list[Timing] = [] + universe = set(range(2**width)) + for density in densities: + on_set = {m for m in universe if rng.random() < density} + if not on_set: + continue + primes_s, chart_s, cover_s, primes = measure_minimisation(on_set, set(), width) + results.append( + Timing( + label=f"density={density:g}", + cases=0, + conditions=width, + positives=len(on_set), + remainders=0, + primes=primes, + truth_table_s=0.0, + primes_s=primes_s, + chart_s=chart_s, + cover_s=cover_s, + ) + ) + return results + + +def sweep_remainders(rng: Random, width: int, shares: tuple[float, ...]) -> list[Timing]: + """Vary how many configurations are remainders, at a fixed on-set.""" + results: list[Timing] = [] + universe = set(range(2**width)) + on_set = {m for m in universe if rng.random() < 0.15} + if not on_set: + return results + rest = sorted(universe - on_set) + for share in shares: + count = int(len(rest) * share) + dont_cares = set(rest[:count]) + primes_s, chart_s, cover_s, primes = measure_minimisation(on_set, dont_cares, width) + results.append( + Timing( + label=f"remainders={share:g}", + cases=0, + conditions=width, + positives=len(on_set), + remainders=len(dont_cares), + primes=primes, + truth_table_s=0.0, + primes_s=primes_s, + chart_s=chart_s, + cover_s=cover_s, + ) + ) + return results + + +def render(results: list[Timing], *, markdown: bool) -> str: + """Format the timings as a plain or Markdown table.""" + header = ( + "label", + "cases", + "cond", + "pos", + "rem", + "primes", + "table_s", + "primes_s", + "chart_s", + "cover_s", + "total_s", + "dominant", + ) + rows = [ + ( + item.label, + str(item.cases), + str(item.conditions), + str(item.positives), + str(item.remainders), + str(item.primes), + f"{item.truth_table_s:.4f}", + f"{item.primes_s:.4f}", + f"{item.chart_s:.4f}", + f"{item.cover_s:.4f}", + f"{item.total_s:.4f}", + item.dominant, + ) + for item in results + ] + if markdown: + lines = ["| " + " | ".join(header) + " |", "| " + " | ".join("---" for _ in header) + " |"] + lines.extend("| " + " | ".join(row) + " |" for row in rows) + return "\n".join(lines) + + widths = [max(len(header[i]), *(len(row[i]) for row in rows)) for i in range(len(header))] + lines = [" ".join(name.rjust(widths[i]) for i, name in enumerate(header))] + lines.append(" ".join("-" * width for width in widths)) + lines.extend(" ".join(cell.rjust(widths[i]) for i, cell in enumerate(row)) for row in rows) + return "\n".join(lines) + + +def main() -> None: + """Run the sweeps and print the results.""" + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument("--max-width", type=int, default=8) + parser.add_argument("--cases", type=int, default=40) + parser.add_argument("--seed", type=int, default=20260811) + parser.add_argument("--markdown", action="store_true") + args = parser.parse_args() + + # Each sweep gets its own generator. Sharing one would make every table + # depend on how much randomness the tables before it happened to consume, + # so --max-width would silently change the numbers further down. + sections = [ + ("Conditions", sweep_conditions(Random(args.seed), args.max_width, args.cases)), + ("Cases", sweep_cases(Random(args.seed + 1), 6, (10, 25, 50, 100, 250))), + ("Sufficient share", sweep_density(Random(args.seed + 2), 7, (0.1, 0.25, 0.5, 0.75))), + ("Remainder share", sweep_remainders(Random(args.seed + 3), 7, (0.0, 0.25, 0.5, 1.0))), + ] + for title, results in sections: + if not results: + continue + print(f"\n### {title}\n" if args.markdown else f"\n=== {title} ===") + print(render(results, markdown=args.markdown)) + + +if __name__ == "__main__": + main() diff --git a/docs/guide/minimisation.md b/docs/guide/minimisation.md index 91c2e21..716d4aa 100644 --- a/docs/guide/minimisation.md +++ b/docs/guide/minimisation.md @@ -67,30 +67,35 @@ implementation escapes that. What matters in practice is not the number of conditions alone but the *shape* of the chart: how many prime implicants there are, and how much they overlap. -The figures below are indicative measurements on ordinary hardware, for a -typical QCA design: 40 observed cases, about a third of the observed rows -sufficient, and every remaining row treated as a logical remainder — that is, -the **parsimonious** solution, which is the more expensive of the two standard -families because it hands the solver a large don't-care set. +The figures below come from `benchmarks/profile_phases.py` on ordinary hardware, +for a typical QCA design: 40 observed cases, an outcome that genuinely depends on +the conditions, and every unobserved row treated as a logical remainder — that +is, the **parsimonious** solution, which is the more expensive of the two +standard families because it hands the solver a large don't-care set. | Conditions | Truth-table rows | Remainders | Parsimonious solution | | --- | --- | --- | --- | -| 6 | 64 | 24 | ~0.004 s | -| 7 | 128 | 88 | ~0.05 s | -| 8 | 256 | 216 | ~0.6 s | -| 9 | 512 | 472 | ~6 s | -| 10 | 1 024 | 984 | ~36 s | +| 6 | 64 | 36 | 0.004 s | +| 7 | 128 | 95 | 0.009 s | +| 8 | 256 | 218 | 0.034 s | +| 9 | 512 | 473 | 0.113 s | +| 10 | 1 024 | 984 | 0.448 s | -Roughly an order of magnitude per additional condition. Conservative solutions -are considerably cheaper at the same width, since they use no don't-cares. +Roughly a trebling per additional condition in this regime, where the cost is +carried by prime generation and chart construction rather than by the exponential +search. Conservative solutions are cheaper still at the same width, since they +use no don't-cares. Dense tables — where a large fraction of *all* minterms is sufficient — are the -worst case for the chart solver and degrade sooner. This is the regime -`benchmarks/benchmark_qmc.py` measures: +worst case, because they are the regime where the exact search itself dominates. +At seven conditions, timing the cover phase alone: -```bash -python benchmarks/benchmark_qmc.py --max-width 9 --density 0.3 -``` +| Sufficient share | Positives | Primes | Cover solving | +| --- | --- | --- | --- | +| 0.10 | 11 | 10 | <0.0001 s | +| 0.25 | 33 | 26 | 0.0001 s | +| 0.50 | 65 | 53 | 0.0024 s | +| 0.75 | 91 | 78 | 0.60 s | Do not extrapolate from either table to your own data; the chart shape matters more than the width. If a run does not finish, the practical levers are reducing @@ -101,6 +106,62 @@ This cost is the price of exactness. A faster minimiser (CCubes/eQMC-style) is a roadmap item, but it will be added as an alternative engine rather than by weakening the guarantee of the current one. +## Being told before it gets slow + +Exact minimisation is worst-case exponential, and that cost is not negotiable +here: no heuristic is substituted when a problem gets hard, because a silently +approximate answer is worse than a slow one. What *is* negotiable is finding out +in advance. + +Once the primes are generated — which is fast — the chart's shape is known, and +`minimize` warns before entering the exponential phase: + +```python +>>> minimize(large_on_set, width=12) +MinimizationComplexityWarning: Exact minimisation of 400 configurations over 12 +conditions produced 180 prime implicants (72000 chart cells). Solving the chart +exactly is worst-case exponential and may take a long time. The result will +still be exact. To reduce the cost, use fewer conditions, tighten the +consistency cutoff so fewer configurations qualify, or lower max_solutions if +the model is highly ambiguous. +``` + +The run still completes, and still returns an exact answer. Silence it with +`warnings.simplefilter`, or switch the check off with `complexity_guard=False`. + +The trigger is the **number of prime implicants**, not the number of conditions: +the density table above is the evidence, and the threshold sits either side of +the climb between 53 and 78 primes. That climb is driven by how much the primes +*overlap*, which is why this is a warning rather than a prediction — a chart with +many primes and little overlap solves instantly. + +## Where the time goes + +`benchmarks/profile_phases.py` times each phase separately across four +dimensions — cases, conditions, sufficient configurations and remainders: + +```bash +python benchmarks/profile_phases.py +python benchmarks/profile_phases.py --max-width 9 --markdown +``` + +Two different phases dominate in two different regimes: + +- **Remainder-heavy problems** — the parsimonious case, where most of the + property space is unobserved — are dominated by prime generation. +- **Dense on-sets** are dominated by solving the chart. + +Truth-table construction barely moves with the number of cases — 25× the cases +costs under 2× the time, because the work is proportional to the property space, +not to the sample. Adding cases is cheap; adding *conditions* is not. + +!!! note "Measure before optimising" + Prime generation was originally the bottleneck at eight conditions, taking + 1.4 s. Rewriting it over integer masks rather than tuples of optional bits + brought that to 0.010 s — the same algorithm, a different representation. + The exactness and R-parity tests passed unchanged, which is what made the + rewrite safe to keep. + ## How the search stays tractable Three reductions cut the search space without ever changing the answer: diff --git a/src/setqca/__init__.py b/src/setqca/__init__.py index f615fbc..ce1b5e0 100644 --- a/src/setqca/__init__.py +++ b/src/setqca/__init__.py @@ -57,6 +57,7 @@ BooleanSolution, Implicant, MinimalCover, + MinimizationComplexityWarning, MinimizationResult, PrimeImplicant, PrimeImplicantChart, @@ -99,6 +100,7 @@ "Implication", "Intersection", "MinimalCover", + "MinimizationComplexityWarning", "MinimizationResult", "MultiValueDomain", "MultiValueResult", diff --git a/src/setqca/minimize/__init__.py b/src/setqca/minimize/__init__.py index cd88922..82b6494 100644 --- a/src/setqca/minimize/__init__.py +++ b/src/setqca/minimize/__init__.py @@ -8,17 +8,25 @@ build_chart, minimize_chart, ) +from .complexity import ( + ComplexityEstimate, + MinimizationComplexityWarning, + estimate_complexity, +) from .implicant import Implicant from .qmc import BooleanSolution, exact_minimum_covers, minimize, prime_implicants __all__ = [ "BooleanSolution", + "ComplexityEstimate", "Implicant", "MinimalCover", + "MinimizationComplexityWarning", "MinimizationResult", "PrimeImplicant", "PrimeImplicantChart", "build_chart", + "estimate_complexity", "exact_minimum_covers", "minimize", "minimize_chart", diff --git a/src/setqca/minimize/chart.py b/src/setqca/minimize/chart.py index 12b7f95..e12ea6a 100644 --- a/src/setqca/minimize/chart.py +++ b/src/setqca/minimize/chart.py @@ -239,8 +239,9 @@ def build_chart( *, dont_cares: set[int] | None = None, width: int, + primes: tuple[Implicant, ...] | None = None, ) -> PrimeImplicantChart: - """Generate the prime implicants and assemble the chart. + """Assemble the chart, generating the prime implicants if not supplied. Parameters ---------- @@ -250,6 +251,10 @@ def build_chart( Logical remainders, usable but not required. width : int Number of conditions. + primes : tuple of Implicant, optional + Pre-generated primes. Supply these when they are already available: + generation dominates the runtime on remainder-heavy problems, so + repeating it doubles the cost of building a chart. Returns ------- @@ -257,7 +262,8 @@ def build_chart( The chart, whose primes are ordered by literal count then bit pattern. """ required = set(on_set) - primes = prime_implicants(required, set() if dont_cares is None else set(dont_cares), width) + if primes is None: + primes = prime_implicants(required, set() if dont_cares is None else set(dont_cares), width) entries = tuple( PrimeImplicant( index=index, @@ -302,7 +308,10 @@ def minimize_chart( RuntimeError If some configuration is covered by no prime implicant. """ - chart = build_chart(on_set, dont_cares=dont_cares, width=width) + generated = prime_implicants( + set(on_set), set() if dont_cares is None else set(dont_cares), width + ) + chart = build_chart(on_set, dont_cares=dont_cares, width=width, primes=generated) solutions = exact_minimum_covers( tuple(prime.implicant for prime in chart.primes), set(on_set), diff --git a/src/setqca/minimize/complexity.py b/src/setqca/minimize/complexity.py new file mode 100644 index 0000000..1ee9cc2 --- /dev/null +++ b/src/setqca/minimize/complexity.py @@ -0,0 +1,153 @@ +"""Warning before combinatorial explosion, rather than after. + +Exact Boolean minimisation is worst-case exponential. That cost is the price of +an exact answer and is not negotiable here — no heuristic is substituted when a +problem gets hard, because a solution that is silently approximate is worse than +one that takes a long time. + +What *is* negotiable is being told in advance. This module estimates how hard a +chart looks before the exponential phase begins, so a run that will not finish +says so rather than appearing to hang. + +Where the thresholds come from +------------------------------ + +Measured with ``benchmarks/profile_phases.py`` on tables of seven conditions, +timing the exact cover phase alone: + +=========== ========== =============== +Primes Positives Cover solving +=========== ========== =============== +10 11 <0.0001 s +26 33 0.0001 s +53 65 0.0024 s +78 91 0.60 s +=========== ========== =============== + +The count of prime implicants is the useful predictor: the cost climbs sharply +somewhere past sixty, and the climb is driven by how much the primes overlap +rather than by the number of conditions on its own. The thresholds below sit +either side of that transition. They are a warning, not a prediction — a chart +with many primes and little overlap solves instantly. +""" + +from __future__ import annotations + +import warnings +from dataclasses import dataclass +from typing import Literal + +MODERATE_PRIMES = 32 +HIGH_PRIMES = 64 + +Level = Literal["low", "moderate", "high"] + + +class MinimizationComplexityWarning(UserWarning): + """Raised when a minimisation problem looks likely to be slow. + + The computation still runs, and still returns an exact answer. Silence it + with :func:`warnings.simplefilter`, or turn the check off entirely with + ``complexity_guard=False``. + """ + + +@dataclass(frozen=True, slots=True) +class ComplexityEstimate: + """How hard a prime-implicant chart looks before it is solved.""" + + width: int + required: int + dont_cares: int + primes: int + + @property + def universe(self) -> int: + """Return the size of the property space.""" + size: int = 2**self.width + return size + + @property + def chart_cells(self) -> int: + """Return the size of the chart, primes times rows to cover.""" + return self.primes * self.required + + @property + def level(self) -> Level: + """Return a coarse difficulty band.""" + if self.primes > HIGH_PRIMES: + return "high" + if self.primes > MODERATE_PRIMES: + return "moderate" + return "low" + + @property + def should_warn(self) -> bool: + """Return whether this problem warrants warning the caller.""" + return self.level == "high" + + @property + def message(self) -> str: + """Return a description of the problem's shape and what to do about it.""" + return ( + f"Exact minimisation of {self.required} configurations over " + f"{self.width} conditions produced {self.primes} prime implicants " + f"({self.chart_cells} chart cells). Solving the chart exactly is " + "worst-case exponential and may take a long time. The result will " + "still be exact. To reduce the cost, use fewer conditions, tighten " + "the consistency cutoff so fewer configurations qualify, or lower " + "max_solutions if the model is highly ambiguous." + ) + + def __str__(self) -> str: + return ( + f"width={self.width}, required={self.required}, " + f"dont_cares={self.dont_cares}, primes={self.primes}, " + f"level={self.level}" + ) + + +def estimate_complexity( + *, width: int, required: int, dont_cares: int, primes: int +) -> ComplexityEstimate: + """Describe how hard a chart looks, from counts alone. + + Parameters + ---------- + width : int + Number of conditions. + required : int + Configurations that must be covered. + dont_cares : int + Configurations usable but not required. + primes : int + Number of prime implicants generated. + + Returns + ------- + ComplexityEstimate + The counts, a difficulty band, and an explanatory message. + """ + return ComplexityEstimate(width=width, required=required, dont_cares=dont_cares, primes=primes) + + +def warn_if_complex(estimate: ComplexityEstimate, *, stacklevel: int = 3) -> bool: + """Emit :class:`MinimizationComplexityWarning` when the chart looks hard. + + Parameters + ---------- + estimate : ComplexityEstimate + The estimate to judge. + stacklevel : int, default 3 + Passed through to :func:`warnings.warn` so the warning points at the + caller's code rather than at this module. + + Returns + ------- + bool + Whether a warning was emitted. + """ + if not estimate.should_warn: + return False + warnings.warn(estimate.message, MinimizationComplexityWarning, stacklevel=stacklevel) + return True diff --git a/src/setqca/minimize/qmc.py b/src/setqca/minimize/qmc.py index f319b76..8124b00 100644 --- a/src/setqca/minimize/qmc.py +++ b/src/setqca/minimize/qmc.py @@ -6,7 +6,8 @@ from collections.abc import Sequence from dataclasses import dataclass -from .implicant import Implicant, minterm_to_implicant +from .complexity import estimate_complexity, warn_if_complex +from .implicant import Implicant @dataclass(frozen=True, slots=True) @@ -26,47 +27,82 @@ def as_expression(self, conditions: tuple[str, ...]) -> str: def prime_implicants(on_set: set[int], dont_cares: set[int], width: int) -> tuple[Implicant, ...]: - """Generate all prime implicants exactly using classical QMC.""" + """Generate all prime implicants exactly using classical QMC. + + Internally a cube is a pair of integers — a mask of the fixed positions and + the values at those positions — rather than a tuple of optional bits. + Benchmarking showed this phase dominating on remainder-heavy problems, + which is exactly the parsimonious case, and integer masks make combining + a few machine operations instead of a tuple walk. The algorithm is + unchanged; only the representation is. + + Parameters + ---------- + on_set : set of int + Minterms that must be covered. + dont_cares : set of int + Minterms usable but not required. + width : int + Number of conditions. + + Returns + ------- + tuple of Implicant + Primes covering at least one required minterm, ordered by literal + count then bit pattern. + + Raises + ------ + ValueError + If the two sets overlap. + """ if on_set & dont_cares: raise ValueError("on_set and dont_cares must be disjoint.") universe = on_set | dont_cares if not universe: return () - current = {minterm_to_implicant(value, width) for value in universe} - primes: set[Implicant] = set() + full_mask = (1 << width) - 1 + # A cube is (mask, value): `mask` marks the fixed positions, `value` holds + # the bits there. A minterm m is covered when `m & mask == value`. + current: set[tuple[int, int]] = {(full_mask, minterm) for minterm in universe} + primes: set[tuple[int, int]] = set() while current: - grouped: dict[int, list[Implicant]] = defaultdict(list) - for implicant in current: - grouped[sum(bit == 1 for bit in implicant.pattern)].append(implicant) - - used: set[Implicant] = set() - next_round: dict[tuple[int | None, ...], Implicant] = {} - for ones in sorted(grouped): - for left in grouped[ones]: - for right in grouped.get(ones + 1, []): - combined = left.combine(right) - if combined is None: + grouped: dict[tuple[int, int], list[int]] = defaultdict(list) + for mask, value in current: + grouped[(mask, value.bit_count())].append(value) + + used: set[tuple[int, int]] = set() + next_round: set[tuple[int, int]] = set() + for (mask, ones), values in grouped.items(): + neighbours = grouped.get((mask, ones + 1)) + if not neighbours: + continue + for left in values: + for right in neighbours: + difference = left ^ right + # Exactly one differing fixed position, and both cubes fix + # it, is the classical combination rule. + if difference & (difference - 1) or not difference & mask: continue - used.add(left) - used.add(right) - previous = next_round.get(combined.pattern) - if previous is None: - next_round[combined.pattern] = combined - else: - next_round[combined.pattern] = Implicant( - combined.pattern, previous.origins | combined.origins - ) - - primes.update(item for item in current if item not in used) - current = set(next_round.values()) + used.add((mask, left)) + used.add((mask, right)) + reduced = mask & ~difference + next_round.add((reduced, left & reduced)) + + primes.update(cube for cube in current if cube not in used) + current = next_round + + def covers(cube: tuple[int, int], minterm: int) -> bool: + mask, value = cube + return minterm & mask == value # A prime built only from don't-cares cannot cover any required minterm. - useful = [item for item in primes if any(item.covers(m) for m in on_set)] + useful = [cube for cube in primes if any(covers(cube, m) for m in on_set)] return tuple( sorted( - useful, + (_to_implicant(cube, width, on_set | dont_cares) for cube in useful), key=lambda item: ( item.literals, tuple(2 if bit is None else bit for bit in item.pattern), @@ -75,6 +111,16 @@ def prime_implicants(on_set: set[int], dont_cares: set[int], width: int) -> tupl ) +def _to_implicant(cube: tuple[int, int], width: int, universe: set[int]) -> Implicant: + """Convert an internal (mask, value) cube to the public implicant type.""" + mask, value = cube + pattern = tuple( + None if not mask >> shift & 1 else value >> shift & 1 for shift in reversed(range(width)) + ) + origins = frozenset(minterm for minterm in universe if minterm & mask == value) + return Implicant(pattern, origins) + + def exact_minimum_covers( primes: tuple[Implicant, ...], on_set: set[int], @@ -230,8 +276,40 @@ def minimize( dont_cares: set[int] | None = None, width: int, max_solutions: int = 256, + complexity_guard: bool = True, ) -> tuple[BooleanSolution, ...]: - """Return all exact minimum Boolean covers for the specified truth table.""" + """Return all exact minimum Boolean covers for the specified truth table. + + Parameters + ---------- + on_set : set of int + Minterms that must be covered. + dont_cares : set of int, optional + Logical remainders, usable but not required. + width : int + Number of conditions. + max_solutions : int, default 256 + Upper bound on the number of tied minimum covers returned. + complexity_guard : bool, default True + Warn with :class:`~setqca.minimize.MinimizationComplexityWarning` when + the chart looks likely to be slow. The result is exact either way; the + warning arrives before the expensive phase rather than after it. + + Returns + ------- + tuple of BooleanSolution + Every cover of provably minimal cost. + """ dc = set() if dont_cares is None else set(dont_cares) - primes = prime_implicants(set(on_set), dc, width) - return exact_minimum_covers(primes, set(on_set), max_solutions=max_solutions) + required = set(on_set) + primes = prime_implicants(required, dc, width) + if complexity_guard: + warn_if_complex( + estimate_complexity( + width=width, + required=len(required), + dont_cares=len(dc), + primes=len(primes), + ) + ) + return exact_minimum_covers(primes, required, max_solutions=max_solutions) diff --git a/tests/test_complexity.py b/tests/test_complexity.py new file mode 100644 index 0000000..3c21f43 --- /dev/null +++ b/tests/test_complexity.py @@ -0,0 +1,115 @@ +"""Tests for the minimisation complexity guard.""" + +from __future__ import annotations + +import warnings + +import pytest + +from setqca.minimize import complexity, minimize +from setqca.minimize.complexity import ( + HIGH_PRIMES, + MODERATE_PRIMES, + ComplexityEstimate, + MinimizationComplexityWarning, + estimate_complexity, + warn_if_complex, +) + + +def _estimate(primes: int) -> ComplexityEstimate: + return estimate_complexity(width=8, required=40, dont_cares=10, primes=primes) + + +class TestEstimate: + def test_counts_are_reported_back(self) -> None: + estimate = estimate_complexity(width=6, required=12, dont_cares=30, primes=9) + assert estimate.width == 6 + assert estimate.required == 12 + assert estimate.dont_cares == 30 + assert estimate.primes == 9 + assert estimate.universe == 64 + assert estimate.chart_cells == 9 * 12 + + @pytest.mark.parametrize( + ("primes", "level"), + [ + (1, "low"), + (MODERATE_PRIMES, "low"), + (MODERATE_PRIMES + 1, "moderate"), + (HIGH_PRIMES, "moderate"), + (HIGH_PRIMES + 1, "high"), + ], + ) + def test_the_band_follows_the_prime_count(self, primes: int, level: str) -> None: + assert _estimate(primes).level == level + + def test_only_the_high_band_warrants_a_warning(self) -> None: + assert not _estimate(MODERATE_PRIMES).should_warn + assert not _estimate(HIGH_PRIMES).should_warn + assert _estimate(HIGH_PRIMES + 1).should_warn + + def test_the_message_says_the_result_is_still_exact(self) -> None: + """The guard must not read as a threat to correctness.""" + message = _estimate(HIGH_PRIMES + 1).message + assert "still be exact" in message + assert "fewer conditions" in message + + def test_the_estimate_renders_readably(self) -> None: + assert "level=high" in str(_estimate(HIGH_PRIMES + 1)) + + +class TestWarning: + def test_an_easy_problem_is_silent(self) -> None: + with warnings.catch_warnings(): + warnings.simplefilter("error") + assert warn_if_complex(_estimate(4)) is False + + def test_a_hard_problem_warns(self) -> None: + with pytest.warns(MinimizationComplexityWarning, match="worst-case exponential"): + assert warn_if_complex(_estimate(HIGH_PRIMES + 1)) is True + + def test_the_warning_is_a_user_warning(self) -> None: + """So that ordinary warning filters reach it.""" + assert issubclass(MinimizationComplexityWarning, UserWarning) + + +class TestIntegration: + def test_ordinary_problems_do_not_warn(self) -> None: + """The suite runs with warnings as errors, so this is enforced throughout.""" + with warnings.catch_warnings(): + warnings.simplefilter("error") + minimize({6, 7}, dont_cares={4, 5}, width=3) + + def test_the_guard_fires_before_the_expensive_phase( + self, monkeypatch: pytest.MonkeyPatch + ) -> None: + """Lowering the threshold proves the wiring without paying the cost. + + A chart genuinely large enough to trip the real threshold is, by + construction, one that takes a long time to solve — so the test lowers + the bar instead of building such a chart. + """ + monkeypatch.setattr(complexity, "HIGH_PRIMES", 0) + with pytest.warns(MinimizationComplexityWarning, match="still be exact"): + minimize({6, 7}, width=3) + + def test_the_guard_can_be_switched_off(self, monkeypatch: pytest.MonkeyPatch) -> None: + monkeypatch.setattr(complexity, "HIGH_PRIMES", 0) + with warnings.catch_warnings(): + warnings.simplefilter("error") + minimize({6, 7}, width=3, complexity_guard=False) + + def test_warning_does_not_change_the_answer(self, monkeypatch: pytest.MonkeyPatch) -> None: + quiet = minimize({6, 7}, dont_cares={4, 5}, width=3, complexity_guard=False) + monkeypatch.setattr(complexity, "HIGH_PRIMES", 0) + with pytest.warns(MinimizationComplexityWarning): + loud = minimize({6, 7}, dont_cares={4, 5}, width=3) + assert [s.implicants for s in loud] == [s.implicants for s in quiet] + + def test_the_warning_is_catchable_from_the_top_level_namespace(self) -> None: + """Callers need the class to filter on, so it is exported alongside the API.""" + import setqca + + assert setqca.MinimizationComplexityWarning is MinimizationComplexityWarning + assert "MinimizationComplexityWarning" in setqca.__all__ diff --git a/tests/test_errors.py b/tests/test_errors.py index 1fa51a6..a9326d6 100644 --- a/tests/test_errors.py +++ b/tests/test_errors.py @@ -67,6 +67,19 @@ def test_implicants_of_different_widths_do_not_combine(self) -> None: def test_implicants_differing_in_more_than_one_literal_do_not_combine(self) -> None: assert minterm_to_implicant(0, 3).combine(minterm_to_implicant(3, 3)) is None + def test_adjacent_implicants_combine_into_a_wider_cube(self) -> None: + """`Implicant.combine` is public API, so it is tested on its own merits. + + Prime generation now works on integer masks internally and no longer + calls this, but the method remains part of the type's interface. + """ + combined = minterm_to_implicant(6, 3).combine(minterm_to_implicant(7, 3)) + assert combined is not None + assert combined.pattern == (1, 1, None) + assert combined.origins == frozenset({6, 7}) + assert combined.literals == 2 + assert combined.as_expression(("A", "B", "C")) == "A*B" + def test_a_dont_care_literal_blocks_combination(self) -> None: left = Implicant((1, None, 0), frozenset({4})) right = Implicant((1, 1, 1), frozenset({7}))