Constraint Satisfaction Problems

A Constraint Satisfaction Problem is a triple $(X, D, C)$ where:

• $X = {X_1, X_2, \ldots, X_n}$ is a finite set of variables.
• $D = {D_1, D_2, \ldots, D_n}$ is a set of domains, where $D_i$ is the finite set of possible values for variable $X_i$.
• $C = {C_1, C_2, \ldots, C_m}$ is a set of constraints, where each constraint $C_k$ is defined over some subset (scope) of variables and specifies the legal combinations of values for those variables.

A complete assignment maps every variable to a value from its domain. A complete assignment is a solution if it satisfies all constraints. A partial assignment assigns values to some subset of variables. We say a partial assignment is consistent (or locally consistent) if it violates no constraint involving only already-assigned variables.

Example — Map Coloring (Australia): Variables: WA, NT, SA, Q, NSW, V, T (the states/territories) Domain for each: {red, green, blue} Constraints: WA≠NT, WA≠SA, NT≠SA, NT≠Q, SA≠Q, SA≠NSW, SA≠V, Q≠NSW, NSW≠V (adjacent states must differ)

A solution is any assignment of colors that satisfies all adjacency constraints.

Types of constraints by arity:

• Unary: involves one variable (e.g., $X_1 \neq \text{red}$)
• Binary: involves two variables (e.g., $X_1 \neq X_2$) — the most common type
• Higher-order (global): involves three or more variables (e.g., AllDifferent(X_1, X_2, X_3, X_4) in Sudoku)

Global constraints like AllDifferent are enormously useful because efficient special-purpose propagators can be written for them — far better than decomposing them into binary constraints.

# A simple CSP representation in Python
from collections import defaultdict

class CSP:
    def __init__(self):
        self.variables = []         # list of variable names
        self.domains = {}           # var -> list of values
        self.constraints = defaultdict(list)  # var -> list of (neighbor, predicate)

    def add_variable(self, var, domain):
        self.variables.append(var)
        self.domains[var] = list(domain)

    def add_binary_constraint(self, var1, var2, predicate):
        """predicate(val1, val2) -> bool"""
        self.constraints[var1].append((var2, predicate))
        self.constraints[var2].append((var1, lambda v2, v1: predicate(v1, v2)))

# Australia map-coloring
csp = CSP()
for state in ['WA','NT','SA','Q','NSW','V','T']:
    csp.add_variable(state, ['red','green','blue'])
for u, v in [('WA','NT'),('WA','SA'),('NT','SA'),('NT','Q'),
             ('SA','Q'),('SA','NSW'),('SA','V'),('Q','NSW'),('NSW','V')]:
    csp.add_binary_constraint(u, v, lambda a, b: a != b)

For CSPs involving only binary constraints, the constraint graph is the most illuminating data structure: nodes are variables, and an edge connects two variables if and only if they share a constraint.

The topology of the constraint graph is not just a visualization tool — it directly determines the hardness of solving the CSP:

1. Trees (acyclic constraint graphs): CSPs whose constraint graph is a tree can be solved in $O(nd^2)$ time — polynomial — using a simple directed arc-consistency pass followed by backtrack-free assignment. This is a remarkable result: just by structuring the problem as a tree, exponential search becomes polynomial.

2. Sparse graphs: problems with low treewidth $w$ can be solved in $O(d^w \cdot n)$ time using tree decomposition. Treewidth measures how tree-like the graph is.

3. Dense graphs: in the worst case (complete graph), $O(d^n)$ — genuine exponential.

For CSPs with higher-order constraints, the analogous structure is the constraint hypergraph (or the dual graph, where constraint nodes and variable nodes both appear). In practice, Sudoku's AllDifferent constraints each span 9 variables; modeling them as binary constraints would produce $\binom{9}{2}=36$ binary constraints per row/column/box, obscuring structure that specialized propagators can exploit.

Why this matters: When you model a new problem as a CSP, one of the first questions to ask is: can I choose a variable representation that makes the constraint graph sparse or tree-like? This is not just an academic exercise — industrial constraint solvers routinely exploit tree decomposition to tackle problems with thousands of variables.

# Tree-structured CSP: O(n*d^2) exact solver
# Topological sort + directed AC + backtrack-free assignment

def tree_csp_solver(csp, root):
    """Assumes csp.constraint_graph is a tree."""
    order = topological_sort(csp, root)   # BFS/DFS order from root

    # Phase 1: directed arc consistency, leaves -> root
    for Xi in reversed(order):
        parent = get_parent(order, Xi)
        if parent is not None:
            revise(csp, parent, Xi)  # remove values from parent's domain
            if not csp.domains[parent]:
                return None  # no solution

    # Phase 2: assign root-to-leaves, any consistent value
    assignment = {}
    for Xi in order:
        for val in csp.domains[Xi]:
            # Check consistency with already-assigned parent
            parent = get_parent(order, Xi)
            if parent is None or is_consistent(Xi, val, parent,
                                               assignment[parent], csp):
                assignment[Xi] = val
                break
    return assignment

The conceptually simplest CSP solver is backtracking search — a depth-first search through the space of partial assignments that abandons (backtracks from) any partial assignment the moment it is found to violate a constraint.

The key insight over plain exhaustive search is incrementality: we commit to one variable's value at a time, and check constraints before completing the assignment. If $X_i = v$ already violates a constraint with a previously assigned variable $X_j$, we prune the entire subtree of assignments that extend $X_i = v$ — without ever exploring it.

Worst case is still $O(d^n)$, but in practice backtracking is dramatically better than the $O(n! \cdot d^n)$ of enumerate-then-check.

The algorithm is elegantly recursive:

function BACKTRACK(assignment, csp):
    if assignment is complete: return assignment
    var = SELECT-UNASSIGNED-VARIABLE(csp, assignment)
    for each value in ORDER-DOMAIN-VALUES(var, assignment, csp):
        if value is consistent with assignment:
            assignment[var] = value
            result = BACKTRACK(assignment, csp)
            if result != failure: return result
            remove var from assignment
    return failure

The two hooks — SELECT-UNASSIGNED-VARIABLE and ORDER-DOMAIN-VALUES — are where all the intelligence lives. Naive versions choose arbitrarily; smart heuristics (MRV, LCV) can reduce the search tree by orders of magnitude. This is the modularity that makes backtracking extensible.

Analogy: Imagine filling in a Sudoku by randomly picking an empty cell and trying values one by one, erasing when you hit a contradiction. Backtracking is exactly that, formalized. What distinguishes expert human solvers from this naive algorithm is that humans propagate consequences of each choice — they eliminate candidates from related cells immediately. That is the insight behind constraint propagation.

def backtracking_search(csp):
    return backtrack({}, csp)

def backtrack(assignment, csp):
    if len(assignment) == len(csp.variables):
        return assignment  # complete assignment

    var = select_unassigned_variable(csp, assignment)  # naive: first unassigned
    for value in order_domain_values(var, assignment, csp):  # naive: domain order
        if is_consistent(var, value, assignment, csp):
            assignment[var] = value
            result = backtrack(assignment, csp)
            if result is not None:
                return result
            del assignment[var]  # backtrack
    return None  # failure

def is_consistent(var, value, assignment, csp):
    for neighbor, predicate in csp.constraints[var]:
        if neighbor in assignment:
            if not predicate(value, assignment[neighbor]):
                return False
    return True

def select_unassigned_variable(csp, assignment):
    for v in csp.variables:
        if v not in assignment:
            return v

def order_domain_values(var, assignment, csp):
    return csp.domains[var]

Forward checking is the simplest form of constraint propagation integrated into backtracking: whenever we assign a value to variable $X_i$, we immediately look ahead at each unassigned neighbor $X_j$ in the constraint graph and remove from $D_j$ any value that would be inconsistent with the current assignment of $X_i$.

If any domain is wiped out (becomes empty), we detect failure immediately — before we waste time assigning all other variables — and backtrack.

Example (map coloring):

• Assign WA = red. Forward checking removes red from NT's and SA's domains.
• Assign Q = green. Forward checking removes green from NT's, NSW's domains.
• Now SA = {blue} (only one value left) because red and green were eliminated.
• Assign NT = blue (only consistent remaining value). Forward checking removes blue from SA — SA's domain is now empty! We detect this and backtrack immediately, saving the effort of assigning Q, NSW, and V.

Key insight: Forward checking maintains a weaker form of consistency (arc consistency from assigned to unassigned variables) but does so cheaply at each assignment step. It is the most common enhancement to basic backtracking because it is simple to implement and very effective in practice.

Limitation: Forward checking only looks one step ahead. It does not propagate implications between unassigned variables. Arc consistency (AC-3) fixes this by propagating throughout the constraint network.

def backtrack_fc(assignment, domains, csp):
    """Backtracking with forward checking."""
    if len(assignment) == len(csp.variables):
        return assignment

    var = select_unassigned_variable(csp, assignment, domains)
    for value in order_domain_values(var, assignment, csp):
        if is_consistent(var, value, assignment, csp):
            assignment[var] = value
            # Forward check: prune neighbors' domains
            pruned = {}  # what we removed, for restoration
            failure = False
            for neighbor, predicate in csp.constraints[var]:
                if neighbor not in assignment:
                    pruned[neighbor] = []
                    for nval in list(domains[neighbor]):
                        if not predicate(value, nval):
                            domains[neighbor].remove(nval)
                            pruned[neighbor].append(nval)
                    if not domains[neighbor]:  # domain wipeout!
                        failure = True
                        break
            if not failure:
                result = backtrack_fc(assignment, domains, csp)
                if result is not None:
                    return result
            # Restore pruned values
            for neighbor, vals in pruned.items():
                domains[neighbor].extend(vals)
            del assignment[var]
    return None

Arc consistency is the gold standard of local consistency for binary CSPs. An arc $(X_i, X_j)$ is arc consistent if for every value $x \in D_i$ there exists at least one value $y \in D_j$ such that the constraint $C(X_i, X_j)$ is satisfied. We write this: $\forall x \in D_i, \exists y \in D_j: C(X_i, X_j)(x, y)$.

A CSP is arc consistent (AC) if every arc in both directions is arc consistent.

The canonical algorithm to enforce arc consistency is AC-3:

function AC-3(csp) -> bool:
    queue = all arcs (Xi, Xj) in csp
    while queue is not empty:
        (Xi, Xj) = DEQUEUE(queue)
        if REVISE(csp, Xi, Xj):        # domain of Xi was reduced
            if Di is empty: return false   # contradiction
            for each Xk neighbor of Xi (Xk != Xj):
                ENQUEUE(queue, (Xk, Xi))   # Xk may now lose support
    return true

function REVISE(csp, Xi, Xj) -> bool:
    revised = false
    for each x in Di:
        if no y in Dj satisfies constraint(Xi=x, Xj=y):
            remove x from Di
            revised = true
    return revised

Complexity: AC-3 runs in $O(ed^3)$ time where $e$ is the number of arcs and $d$ is the maximum domain size. Each arc $(X_i, X_j)$ can be enqueued at most $d$ times (once per value removed from $D_j$). Each REVISE call costs $O(d^2)$. Total: $O(ed^3)$.

Why re-enqueue neighbors? When we remove a value from $D_i$, some values in $D_k$ (a neighbor of $X_i$) that previously had support (a consistent partner in $D_i$) may no longer have support. We must recheck those arcs.

AC-3 vs Forward Checking: Forward checking only enforces arcs from the just-assigned variable to unassigned neighbors. AC-3 propagates through the entire network until quiescence. It is strictly stronger — it detects more failures earlier — at the cost of more work per node in the search tree.

Important caveat: Arc consistency is not a complete algorithm. Even after enforcing full arc consistency, the resulting CSP may have no solution (arc consistency cannot detect all inconsistencies). It is a preprocessing and propagation technique used alongside search.

from collections import deque

def ac3(csp):
    """Enforce arc consistency. Returns False if inconsistency detected."""
    # Initialize queue with all arcs in both directions
    queue = deque()
    for var in csp.variables:
        for neighbor, _ in csp.constraints[var]:
            queue.append((var, neighbor))

    while queue:
        xi, xj = queue.popleft()
        if revise(csp, xi, xj):
            if not csp.domains[xi]:  # domain wipeout
                return False
            # Re-enqueue all arcs pointing INTO xi (except from xj)
            for xk, _ in csp.constraints[xi]:
                if xk != xj:
                    queue.append((xk, xi))
    return True

def revise(csp, xi, xj):
    """Remove values from Di that have no support in Dj."""
    revised = False
    for x in list(csp.domains[xi]):
        # x has support if there exists y in Dj satisfying constraint
        has_support = any(
            pred(x, y)
            for (neighbor, pred) in csp.constraints[xi]
            if neighbor == xj
            for y in csp.domains[xj]
        )
        if not has_support:
            csp.domains[xi].remove(x)
            revised = True
    return revised

The order in which we choose variables to assign in backtracking search has a profound effect on performance. Two complementary heuristics are the cornerstone of modern CSP solvers:

Minimum Remaining Values (MRV) — "Fail-First"

Also called the most-constrained variable heuristic: choose the variable with the fewest legal values remaining in its domain.

Why this is brilliant: If variable $X_i$ has only 1 legal value left, we have no choice — we must assign it next (or detect failure immediately). Assigning the most constrained variable first pushes failures as high as possible in the search tree, pruning the largest subtrees.

Intuitively: if you have a Sudoku cell that can only be 7, fill it in now — don't waste branching on cells that have 9 options.

Degree Heuristic (tiebreaker)

When multiple variables are tied for MRV, break ties by choosing the variable involved in the most constraints on other unassigned variables. This variable is most likely to restrict future choices, so assigning it first catches inconsistencies earlier.

MRV + Degree together form a powerful pair: MRV minimizes the branching factor at each level of the search tree; Degree breaks ties in the most informed way.

The key theoretical point: These heuristics exploit the structure of the specific problem instance (the current domain sizes and constraint graph), not just the problem class. This is why they are so much more effective than static orderings.

Empirical impact: On random CSPs near the phase transition (where problems are hardest), MRV alone can reduce search by 1–3 orders of magnitude compared to static variable ordering.

def mrv_select(csp, assignment, domains):
    """Minimum Remaining Values: pick most constrained unassigned variable."""
    unassigned = [v for v in csp.variables if v not in assignment]
    # Primary: fewest remaining values; Secondary: highest degree
    return min(
        unassigned,
        key=lambda v: (
            len(domains[v]),
            -degree(csp, v, assignment)  # negative = higher degree wins
        )
    )

def degree(csp, var, assignment):
    """Number of constraints on unassigned neighbors."""
    return sum(
        1 for neighbor, _ in csp.constraints[var]
        if neighbor not in assignment
    )

# In the backtracking search, replace:
#   var = select_unassigned_variable(csp, assignment)
# with:
#   var = mrv_select(csp, assignment, domains)

While MRV chooses which variable to assign next, the Least Constraining Value (LCV) heuristic chooses which value to try first for the chosen variable.

Rule: Order values by the number of values they eliminate from the domains of neighboring unassigned variables. Try the value that eliminates the fewest choices for neighbors first.

Intuition: We want to leave the search as flexible as possible. Assigning a value that wipes out choices for neighbors is likely to lead to failure; assigning a value that leaves neighbors with many options gives the search more room to maneuver.

Contrast with MRV: MRV is a fail-first strategy for variable selection (expose failures early). LCV is a succeed-first strategy for value selection (maximize the chance that the current branch leads to a solution). They are complementary: MRV for variables, LCV for values.

When to use LCV: LCV is most beneficial when we expect to find a solution (rather than enumerate all solutions) and when the problem has many solutions. If we only need one solution, finding it quickly matters; LCV points toward the solution-rich part of the space.

Cost: LCV requires computing the number of eliminations for each candidate value, which takes $O(d^2)$ per value per constraint in the worst case. This overhead is worth it for heavily constrained problems but may not be worth it for easy problems.

Real-world analogy: Imagine assigning jobs to workers. The LCV principle says: try the assignment that leaves the most future scheduling flexibility — don't lock in a choice that strands other jobs without workers.

def lcv_order(var, assignment, domains, csp):
    """Least Constraining Value: order values by how few choices they eliminate."""
    def count_eliminations(value):
        total = 0
        for neighbor, predicate in csp.constraints[var]:
            if neighbor not in assignment:
                for nval in domains[neighbor]:
                    if not predicate(value, nval):
                        total += 1
        return total

    return sorted(domains[var], key=count_eliminations)

# Full combined solver: MRV variable selection + LCV value ordering + AC3 propagation
def smart_backtrack(assignment, domains, csp):
    if len(assignment) == len(csp.variables):
        return assignment

    var = mrv_select(csp, assignment, domains)
    for value in lcv_order(var, assignment, domains, csp):
        if is_consistent(var, value, assignment, csp):
            assignment[var] = value
            saved_domains = {v: list(d) for v, d in domains.items()}
            domains[var] = [value]  # fix this variable's domain

            if ac3(csp):  # propagate; modifies domains in place
                result = smart_backtrack(assignment, domains, csp)
                if result is not None:
                    return result

            # Restore domains on backtrack
            for v in domains:
                domains[v] = saved_domains[v]
            del assignment[var]
    return None

Many practical CSPs have constraints that span more than two variables. Rather than decomposing these into binary constraints (which loses structure), modern solvers represent them as global constraints with dedicated propagators.

The most important global constraint is AllDifferent($X_1, X_2, \ldots, X_k$): all $k$ variables must take distinct values. It appears in:

• Sudoku: each row, column, and box must have all different digits.
• Graph coloring (via reformulation).
• Assignment problems: no two entities get the same resource.

Naive decomposition: $\binom{k}{2}$ binary ≠ constraints. For $k=9$ (Sudoku), that's 36 binary constraints per row/column/box — totaling 324 binary constraints for the full puzzle.

Hall's theorem and AllDifferent propagation: A far stronger propagator uses Hall's theorem: AllDifferent on $k$ variables is infeasible if and only if there exists a subset $S$ of the variables such that the union of their domains satisfies $|\bigcup_{X \in S} D_X| < |S|$ (fewer values available than variables to fill). Detecting this requires checking all $2^k$ subsets in general — but efficient algorithms (based on bipartite matching) can enforce generalized arc consistency for AllDifferent in $O(k^{1.5})$ time.

Example: Variables ${A, B, C}$ with $D_A = D_B = D_C = {1, 2}$. No binary constraint fires (each pair can be satisfied). But $|\bigcup D_i| = 2 < 3 = |{A,B,C}|$ — the global propagator immediately detects failure, which the binary constraints missed!

This illustrates a general principle: global constraints with specialized propagators are strictly stronger than their binary decompositions and often dramatically so.

# AllDifferent via bipartite matching (simplified sketch)
# Full implementation uses Hopcroft-Karp matching

def propagate_alldifferent(variables, domains):
    """
    Enforce generalized arc consistency for AllDifferent.
    Remove any value v from domain[X] if v cannot appear in any
    maximum matching of the bipartite graph (variables x values).
    Returns False if infeasible.
    """
    # Build bipartite graph: variables on left, values on right
    # Edge (Xi, v) if v in domains[Xi]
    # Find maximum matching M
    matching = find_max_bipartite_matching(variables, domains)

    if len(matching) < len(variables):
        return False  # Can't assign all variables distinct values

    # Remove values not in any maximum matching
    # (Uses Tarjan's SCC on residual graph — omitted for brevity)
    for var in variables:
        for val in list(domains[var]):
            if not in_some_max_matching(var, val, matching, variables, domains):
                domains[var].remove(val)
                if not domains[var]:
                    return False
    return True

# Hall's theorem check (naive O(2^k), for illustration only)
def hall_check(variables, domains):
    from itertools import combinations
    k = len(variables)
    for size in range(1, k+1):
        for subset in combinations(variables, size):
            union_vals = set()
            for v in subset:
                union_vals |= set(domains[v])
            if len(union_vals) < size:
                return False  # Hall condition violated
    return True

The most profound insight in CSP theory is that the structure of the constraint graph determines the computational complexity of solving the CSP. Not just practically, but in a precise theoretical sense.

Tree-Structured CSPs: $O(nd^2)$

If the constraint graph is a tree, the CSP can be solved in time $O(nd^2)$ without any backtracking:

1. Root the tree at any variable.
2. Apply directed arc consistency from leaves to root: for each node $X_i$ with parent $X_j$, remove from $D_i$ any value with no consistent value in $D_j$.
3. Assign values root-to-leaves: each variable has at least one consistent value (guaranteed by step 2).

Why does this work? In a tree, there is exactly one path between any two variables. The directed AC pass ensures that every partial assignment from root to a node is extendable — no dead ends.

Nearly-Tree-Structured CSPs: Cutset Conditioning

For graphs with a small cycle cutset — a set $S$ of variables whose removal makes the graph acyclic — we can:

1. Enumerate all consistent assignments to $S$ (cost: $O(d^{|S|})$).
2. For each, fix those variables' domains to the assigned value and solve the remaining tree-structured CSP in $O(nd^2)$. Total: $O(d^{|S|} \cdot nd^2)$. If $|S|$ is small (logarithmic in $n$), this is polynomial.

Tree Decomposition and Treewidth

The general framework is tree decomposition: decompose the constraint graph into clusters (subsets of variables) arranged in a tree, where:

• Every variable appears in at least one cluster.
• Every constraint's scope is contained in some cluster.
• For any variable, the set of clusters containing it forms a connected subtree.

The treewidth $w$ is the size of the largest cluster minus 1. CSPs with treewidth $w$ can be solved in $O(n \cdot d^{w+1})$.

Practical takeaway: When modeling a CSP, look for natural tree-like structure. Many real problems (parsing, scheduling on a chain, hierarchical configuration) have low treewidth.

# Tree CSP solver: O(n * d^2), backtrack-free
import networkx as nx  # for tree utilities

def solve_tree_csp(variables, domains, constraints, root):
    """
    Solve a tree-structured CSP in O(n*d^2) without backtracking.
    constraints: dict mapping (Xi, Xj) -> set of allowed (vi, vj) pairs
    """
    # Build tree and get topological order
    G = nx.Graph()
    G.add_nodes_from(variables)
    for (u, v) in constraints:
        G.add_edge(u, v)
    assert nx.is_tree(G), "Constraint graph must be a tree"

    order = list(nx.bfs_tree(G, root).nodes())
    parent = {}
    for i, node in enumerate(order):
        for pred in G.predecessors(node) if G.is_directed() else []:
            parent[node] = pred
    # Re-compute parents from BFS tree
    bfs_tree = nx.bfs_tree(G, root)
    for node in order:
        preds = list(bfs_tree.predecessors(node))
        parent[node] = preds[0] if preds else None

    # Phase 1: directed AC from leaves to root
    for Xi in reversed(order[1:]):
        Xj = parent[Xi]
        key = (Xj, Xi) if (Xj, Xi) in constraints else (Xi, Xj)
        allowed = constraints[key]
        domains[Xj] = [
            vj for vj in domains[Xj]
            if any((vj, vi) in allowed or (vi, vj) in allowed
                   for vi in domains[Xi])
        ]
        if not domains[Xj]:
            return None  # no solution

    # Phase 2: assign root-to-leaves
    assignment = {}
    for Xi in order:
        Xj = parent.get(Xi)
        for val in domains[Xi]:
            if Xj is None or is_consistent_pair(Xi, val, Xj,
                                                assignment[Xj], constraints):
                assignment[Xi] = val
                break
    return assignment

Arc consistency (AC, or 2-consistency) is the most common level of consistency enforced in practice, but it is part of a hierarchy:

• 1-consistency (node consistency): Every value in $D_i$ satisfies all unary constraints on $X_i$. Trivially enforced as a preprocessing step.
• 2-consistency (arc consistency): For every value of $X_i$, there exists a consistent value of $X_j$ for every constraint between them.
• 3-consistency (path consistency): For every consistent assignment to any two variables $X_i, X_j$ and any third variable $X_k$ constrained with them, there exists a value for $X_k$ consistent with $X_i$ and $X_j$.
• k-consistency: For every consistent assignment to any $k-1$ variables, every variable constrained with them can be consistently extended.
• Strong k-consistency: $j$-consistent for all $j \leq k$.

A CSP with $n$ variables that is strongly $n$-consistent can be solved without backtracking by assigning variables one at a time, always finding a consistent extension. This is the theoretical ideal.

Path Consistency (PC-2): Enforced by the PC-2 algorithm. It adds new binary constraints to ensure path consistency. Since it can introduce new constraints (not just prune domains), it is more complex to implement and harder to integrate with search.

Practical tradeoff: Higher consistency levels detect more failures earlier, but the cost per propagation step grows rapidly ($O(n^2 d^3)$ for path consistency vs. $O(e d^3)$ for arc consistency). The consensus in practice is that arc consistency is the sweet spot: it is cheap enough to enforce at every search node and catches the vast majority of failures that consistency detects.

Singleton Arc Consistency (SAC): A more recent consistency notion that is stronger than arc consistency: $X_i = v$ is SAC-consistent if the CSP obtained by fixing $X_i = v$ is arc-consistent. Useful for very hard combinatorial problems.

# Hierarchy of consistency (conceptual, not full implementation)

# Node consistency: O(nd)
def node_consistency(domains, unary_constraints):
    for var, values in domains.items():
        domains[var] = [v for v in values
                        if all(c(v) for c in unary_constraints.get(var, []))]

# Arc consistency (AC-3): O(e * d^3) -- see previous section

# Path consistency enforces: for all consistent (xi=a, xj=b),
# find xk satisfying both. If none exists, add constraint (xi,xj) != (a,b).
# This can create new constraints! O(n^3 * d^3)
def path_consistency_check(xi, a, xj, b, xk, domains, binary_constraints):
    """
    Check if there is a value c in domains[xk] consistent with xi=a and xj=b.
    If not, (xi=a, xj=b) must be removed from the allowed pairs for constraint(xi, xj).
    """
    for c in domains[xk]:
        ok_ik = (not binary_constraints.get((xi, xk))) or \
                ((a, c) in binary_constraints[(xi, xk)])
        ok_jk = (not binary_constraints.get((xj, xk))) or \
                ((b, c) in binary_constraints[(xj, xk)])
        if ok_ik and ok_jk:
            return True  # support found
    return False  # no support: (xi=a, xj=b) must be prohibited

Standard backtracking is chronological: when a dead end is reached, we undo the most recent assignment and try a different value. But the conflict may have been caused by a variable assigned much earlier — undoing the most recent assignment does nothing to resolve it.

Example: Suppose variables are assigned in order $X_1, X_2, X_3, X_4, X_5$ and a conflict arises at $X_5$. The conflict involves $X_5$ and $X_2$ only — $X_3$ and $X_4$ are irrelevant. Chronological backtracking would waste time exploring all values of $X_4$, then all of $X_3$, before finally backtracking to $X_2$.

Backjumping jumps directly to the conflict set — the set of previously assigned variables that caused the failure — instead of backtracking one step at a time.

Conflict-Directed Backjumping (CBJ) (Prosser, 1993):

• Each variable $X_i$ maintains a conflict set $CS(X_i)$: the set of previously assigned variables that are responsible for the current constraint violations at $X_i$.
• When $X_i$ fails (all values in its domain are inconsistent), we jump back to the deepest variable in $CS(X_i)$.
• Crucially, when we jump from $X_i$ to $X_j$, we merge $CS(X_i)$ into $CS(X_j)$, carrying forward the conflict information.

Conflict-Driven Clause Learning (CDCL) in SAT solvers is the modern descendant of this idea — it learns new constraints from conflicts, preventing the same mistake from being repeated. CDCL is the reason modern SAT solvers can handle millions of variables. The connection: CSPs ↔ SAT is deep (every CSP can be encoded as SAT), and the most powerful CSP solvers borrow ideas from CDCL.

No-good learning: Related idea — record the no-good (a conflicting partial assignment) and add it as a new constraint. This turns experience into pruning for all future branches.

def cbj_backtrack(assignment, conflict_sets, csp):
    """
    Conflict-Directed Backjumping.
    conflict_sets[var] = set of previously-assigned variables causing conflict.
    Returns (assignment, jump_target) where jump_target is None on success
    or a variable to jump back to on failure.
    """
    if len(assignment) == len(csp.variables):
        return assignment, None

    var = select_unassigned_variable(csp, assignment)
    conflict_sets[var] = set()

    for value in csp.domains[var]:
        conflicts = get_conflicts(var, value, assignment, csp)
        if conflicts:  # value is inconsistent
            conflict_sets[var] |= conflicts
        else:
            assignment[var] = value
            result, jump_target = cbj_backtrack(assignment, conflict_sets, csp)
            if result is not None:
                return result, None
            if jump_target is not None and jump_target != var:
                # Not our problem — keep jumping
                del assignment[var]
                conflict_sets[var] |= (conflict_sets.get(jump_target, set())
                                       - {var})
                return None, jump_target
            # jump_target == var: we are responsible, try next value
            del assignment[var]

    # All values failed: jump to deepest variable in conflict set
    if conflict_sets[var]:
        jump_to = max(conflict_sets[var],
                      key=lambda v: csp.variables.index(v))
        return None, jump_to
    return None, None  # top-level failure

All algorithms discussed so far are complete: they will find a solution if one exists, or prove none exists. For very large problems where we only need a solution (not proof of unsatisfiability), local search offers a dramatically faster alternative.

The Min-Conflicts algorithm (Minton et al., 1992) is the canonical local search method for CSPs:

1. Initialize a complete assignment (all variables assigned, possibly violating constraints).
2. Repeat until solved (or max-steps exceeded): a. If the assignment satisfies all constraints, return it. b. Pick a randomly chosen conflicted variable $X_i$ (one that participates in at least one violated constraint). c. Assign $X_i$ the value that minimizes the number of violated constraints (breaking ties randomly).

Performance: Min-conflicts solves the $n$-queens problem for $n = 1,000,000$ in about 50 steps on average. This is extraordinary — backtracking would be completely intractable at this scale.

Why does it work so well? The algorithm exploits the density of solutions in the search space. For many over-constrained problems near the satisfiable regime, solutions are dense enough that a greedy gradient descent in constraint-violation space finds them quickly.

Limitations:

• Can get stuck in local minima (assignments where every variable reassignment increases conflicts).
• Incomplete: cannot prove unsatisfiability.
• Performance degrades badly near the phase transition (where problems switch from typically satisfiable to typically unsatisfiable).

In practice: Local search is used for large scheduling problems, hardware verification, and as a subroutine within more sophisticated solvers.

import random

def min_conflicts(csp, max_steps=10000):
    """
    Local search with min-conflicts heuristic.
    Works best when problems have many solutions.
    """
    # Start with a random complete assignment
    assignment = {var: random.choice(csp.domains[var])
                  for var in csp.variables}

    for step in range(max_steps):
        # Find all conflicted variables
        conflicted = [
            var for var in csp.variables
            if count_conflicts(var, assignment[var], assignment, csp) > 0
        ]
        if not conflicted:
            return assignment  # solution found!

        # Pick a random conflicted variable
        var = random.choice(conflicted)

        # Assign the value that minimizes conflicts
        min_val = min(
            csp.domains[var],
            key=lambda v: count_conflicts(var, v, assignment, csp)
        )
        assignment[var] = min_val

    return None  # failed within max_steps

def count_conflicts(var, value, assignment, csp):
    return sum(
        1 for neighbor, pred in csp.constraints[var]
        if neighbor in assignment and not pred(value, assignment[neighbor])
    )