herculesabqp 0.1.2

A convex box-constrained quadratic programming solver with warm starts and active-set polishing.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82

import numpy as np
import scipy.sparse as sp

import herculesabqp


def make_sparse_ata_problem(m: int = 80, n: int = 40, density: float = 0.12, seed: int = 0):
    rng = np.random.default_rng(seed)

    # Build a sparse rectangular matrix A, then form Q = A^T A + ridge * I.
    a = sp.random(
        m,
        n,
        density=density,
        data_rvs=lambda k: rng.normal(size=k),
        format="csr",
    )
    ridge = 1e-2
    q = (a.T @ a + ridge * sp.eye(n, format="csr")).tocsr()

    # Choose a linear term that creates a mixed active set.
    c = rng.normal(loc=-0.15, scale=0.2, size=n)
    lb = np.zeros(n, dtype=float)
    ub = np.ones(n, dtype=float)
    return q, c, lb, ub


def solve_root_and_child():
    q, c, lb, ub = make_sparse_ata_problem()

    # Prepare once when Q and c stay fixed across many solves.
    prepared = herculesabqp.PreparedSolver(
        q,
        c,
        assume_symmetric=True,
        scaling="hessian_diag",
        lipschitz_method="auto",
    )

    root = prepared.solve(
        lb,
        ub,
        max_iter=2_000,
        tol=1e-7,
    )

    x_root = np.asarray(root["x"], dtype=float)
    print("root objective:", root["objective"])
    print("root iterations:", root["iterations"])
    print("root rel_gap:", root["rel_gap"])
    print("root kkt_inf:", root["kkt_inf"])

    # Mimic a branch-and-bound child by fixing a few variables that are already
    # close to the upper bound in the root solution.
    child_lb = lb.copy()
    child_ub = ub.copy()
    upper_active = np.flatnonzero(x_root >= 0.95)
    branch_count = min(4, upper_active.size)
    fixed = upper_active[:branch_count]
    child_lb[fixed] = 1.0
    child_ub[fixed] = 1.0

    child = prepared.solve(
        child_lb,
        child_ub,
        x0=x_root,
        max_iter=500,
        tol=1e-7,
    )

    x_child = np.asarray(child["x"], dtype=float)
    print()
    print("fixed indices:", fixed.tolist())
    print("child objective:", child["objective"])
    print("child iterations:", child["iterations"])
    print("child rel_gap:", child["rel_gap"])
    print("child kkt_inf:", child["kkt_inf"])
    print("child fixed values:", x_child[fixed].tolist())


if __name__ == "__main__":
    solve_root_and_child()