Cardinality-constrained portfolio selection asks for the lowest-risk portfolio that holds at most K names out of n, meets a return target and respects a minimum and maximum size for each holding. It is a mixed-integer quadratic program, and there are many good ways to find good solutions to it.

Our question is a different one. Once a portfolio is in hand, how much can be said about how far it might be from the best possible one — in a form that someone who trusts nothing but arithmetic can check?

A branch-and-bound solver does report a bound. But a reported bound is a number printed by a program: to check it you run the program again and inherit its floating-point tolerances, its licence and its bugs. We wanted a bound that comes with a witness — small enough to send alongside the portfolio, and simple enough that a reader can verify it with code they wrote themselves.

Here x holds the weights, y records which names are held, ℓ are minimum “buy-in” sizes and u are caps.

A certificate is a claimed lower bound D on the best achievable variance, plus a handful of numbers from which a verifier holding the problem data can confirm that every feasible portfolio has variance at least D. The distance between your portfolio and D is then a guarantee, not an estimate. We asked three things of it:

Verdicts are PROVED, REFUTED, NOT PROVED and MODEL_MISMATCH. A floating-point checker cannot establish a comparison with zero slack, so a failure to prove is reported as exactly that, never as a refutation. Every checker answers every input, malformed ones included, with a verdict rather than a crash.

The bound itself is classical: split the covariance as Q = A + diag(d) with A positive semidefinite and d ≥ 0, treat the diagonal part per asset, and you get the perspective relaxation. Its Lagrangian dual is a closed-form expression in a few multipliers — except for one term, a lower bound on the minimum of zᵀAz − wᵀz, which normally means inverting A and trusting whoever did it.

One line of algebra removes that step:

So any point — an exact optimum of the relaxation, a single-precision GPU iterate, or a guess — yields a valid witness. The accuracy of x̂ decides how tight the bound is, never whether it holds. That is what lets everything upstream of the checker be ordinary numerical software, run as fast and as loosely as we like.

The same lemma makes it cheap to go further and prove optima. Fix some names as held and others as excluded, and the same dual, restricted to the remaining names, bounds every portfolio consistent with that choice. A branch-and-bound search therefore needs no new matrix fact at any node: A ≥ 0 is established once at the root, and each pruned node records only which names it fixed, a tangent point and three scalars. The search writes its whole tree to a file, and a separate tree checker — about 270 lines of Python sharing no code with the search — rebuilds every node, confirms the leaves cover every possible support, and re-evaluates each bound with a-priori error bounds. Its verdict is OPTIMAL only when the certified bound is within 10⁻⁶ (relative) of a feasible portfolio.

The tree search is a measuring instrument, not something a user needs: it is how we know how far the root certificate sits from the truth.

Python (NumPy only) evaluates every inequality with outward error bounds and establishes that A is positive semidefinite by a Cholesky factorisation of a slightly shifted matrix, whose computed factor certifies semidefiniteness once backward error is accounted for. It also runs as a standalone script.

C99, written separately and sharing no code with the Python checker, uses only the C standard library.

Exact: where n ≤ 100, the check is repeated in rational arithmetic on the standard library alone.

The checkers import nothing from the prover side, and a test enforces it. At 2,901 assets, checking takes about 0.9 s in Python and 2.0 s in C. The two bound checkers are each under 1,000 lines, deliberately, so that they can be read.

We fixed the measurements before changing any method, because several early “findings” turned out to be artifacts of measurement (more on that below). The headline metric is the certified gap, (f_P − D) / f_P, against the best validated feasible portfolio f_P.

The problems

A registry of 257 problems: 20 from OR-Library, 24 from Cesarone et al., 24 Frangioni–Gentile problems, and 189 real-equity mandates. The mandates cross three US universes (n = 1,437, 1,464 and 2,901), seven rolling two-year windows and nine mandate shapes — index tracking with K = 10 to 50, concentrated (K = 30), diversified active (K = 100), closet indexing (K = 300) and enhanced index (K = 750). The data are weekly adjusted closes of the 2,901 NASDAQ- and NYSE-listed common stocks with a complete five-year history to September 2026.

Many estimates of the same risk

The same returns are re-estimated several ways with everything else held fixed — a raw sample covariance, statistical factor models, two shrinkage estimators and a published fundamental model — because, as it turned out, the risk model decides how tight a certificate can be. The full list is below.

Both sides of the gap are ours

Every gap is taken against portfolios we produced and re-validated in double precision, row by row. An optimum counts as proved only when our search proves it and the independent checker accepts its tree. The open-source solver SCIP runs alongside as a reference; its portfolios, re-validated, enter the best known, but its “optimal” status is never counted as a proof.

“Portfolio optimisation” covers a lot of different problems, and a certificate is only as general as the models it has actually been run on. Here is the full inventory: three problem types, and for the main one, seven ways of estimating risk across nine mandate shapes.

Three problem types

Seven estimates of risk, for the mean-variance problem

Every equity row uses the same weekly returns of up to 2,901 US stocks; only the risk estimate changes. The library builds every one of them without BLAS, so the same returns give byte-identical models — and the same fingerprints — on any machine.

The published benchmark covariances are the field’s standard test sets: OR-Library’s five market indices (Hang Seng, DAX, FTSE 100, S&P 100, Nikkei 225), six real-return universes from Cesarone et al., and the Frangioni–Gentile problems at n = 200, 300 and 400 with buy-in thresholds. Both checkers prove all 68 certificates, and the 32 with n ≤ 100 are also proved in exact rational arithmetic. On the 44 benchmark problems with n ≤ 225, the copositive extension brings the median gap to 0.20% and the worst to 3.5%. Separately, we built five singular problems with 200 to 400 assets from the equity returns to study that extension (below).

Nine mandate shapes

Each equity risk model is crossed with three universes (1,437, 1,464 and 2,901 stocks), seven rolling two-year windows and nine mandates: index tracking with K = 10, 20, 30, 40 and 50; concentrated (K = 30); diversified active (K = 100); closet indexing (K = 300); and enhanced index (K = 750). The sections below take the models in turn.

On the ten-factor tier, the certified gap is 0.17% in the median and 2.13% at worst over all 189 mandates, and producing and checking the certificate takes 5.2 s in the median. Both checkers prove every certificate.

To split that gap into the bound’s share and the portfolio’s, we proved the optima: 181 of 189, with trees the independent checker accepts (179 within ten minutes each, two more given an hour; the other eight stop with a certified gap of at most 0.10%). Against those optima, the portfolio is within 0.012% in the median — exactly optimal on 40 — and the bound within 0.10%. Most of the certified gap is the bound’s, and it is small.

SCIP, given the same models and five minutes each, reports optimality on 183 of the 189. Where both of us report an optimum (175 problems) we agree, and none of its re-validated portfolios lies below a bound we certified.

The single largest effect on how tight the certificate is turned out to be not the method but the covariance matrix fed to it. The same 189 mandates under four estimates of the same risk:

The obvious explanation — more asset-specific risk on the diagonal, a tighter bound — does not hold: identity and constant-correlation shrinkage put similar shares on the diagonal, yet their gaps at K = 10 differ by a factor of 25. Where gaps exceed 1%, about three quarters of the gap is risk the relaxation removes by spreading its weight over more than K names, which no K-name portfolio can do.

That gives a predictor every certificate computes anyway: the share of the relaxation’s weight outside its K largest positions. We fixed thresholds on 162 problems and tested them on 405 held out. Where that share is below 2%, all 268 held-out problems had a gap under 1%. Its rank correlation with the gap is 0.83, against −0.20 for the diagonal share. The library reports it as result.spread.

We also tested a published, open-source forty-factor fundamental model (a market factor, thirty industries and nine styles) on the 27 mandates of its matching window: 0.056% median, under 1% on 24 of 27. One model on one date is a single data point, not a survey.

With a raw sample covariance from 104 weekly returns, the certificate is loose for small K: below K/n = 0.02 the median gap is 32% and the worst 78%. For K/n ≥ 0.15 it stays within 2.7% everywhere.

Our first thought was that our method was simply weak. So we gave SCIP ten minutes on each of 18 of these problems. It left a median gap of 31% against 34% for our certificate, which takes about a second. Its bound is tighter than ours on all 18 — by a few points in most cases, by up to 9 points at worst. Ours is weaker than a good solver’s here, but not by an order of magnitude: with this estimate of risk, the problems themselves are hard to bound.

The reason is structural. If Qv = 0, every admissible diagonal split must be zero wherever v is non-zero (since 0 ≤ vᵀ(Q − diag(d))v = −Σd_iv_i²). With 104 observations of 1,437 assets the null space has dimension at least 1,334 and, in every singular problem we have, touches every asset — so the split is zero, and the bound falls back to the plain continuous relaxation.

We first wrote that no convex certificate could do better on singular problems. That was too strong. The ceiling covers bounds that hold for every vector — but portfolios are non-negative and sum to one, and a bound only has to hold on those.

The right-hand side is again a semidefinite quadratic plus a diagonal, so the same relaxation and the same tangent witness apply unchanged — but off-diagonal mass in G now buys a positive split on a singular matrix. On singular problems of 200 to 400 assets built from the same returns, it closes 87–94% of the root gap (for example 11.2% → 0.75% at n = 200), with no upward drift across that range.

The catch is finding G, a semidefinite program over an n × n matrix. Our GPU solver tracks the family’s ceiling up to about n = 250; there, on the 44 registry problems of that size, it cuts the worst gap from 23.8% to 3.5% and the median from 1.30% to 0.20%, improving 20 problems and worsening none, at a median 9.7 s extra. Beyond that it stops converging usefully. The large singular problems are, we now think, hard for this certificate because of the solver, not because of convexity.

The tangent witness is specific to quadratic risk. For the other two problem types the route to a certificate is different, and so is what it can promise.

Index tracking: proving the global optimum

Pick at most K of an index’s n constituents, each weight a whole number of lots between a minimum and a maximum, to minimise tracking error against the index. Here relaxation is not weak but hopeless. On instances with 60 to 98 assets, the root relaxation sits 72–97% below the proven optimum, because fewer names genuinely track worse: 98 names track the index about 24 times better than 10. Any bound that drops the cardinality limit inherits that factor.

So the certificate accounts for every support instead. For a fixed set of names, every lot-grid portfolio lies inside a continuous capped simplex, so the continuous minimum over that set bounds its best lot-grid portfolio. If no set can beat the incumbent and the incumbent is feasible, the incumbent is the global optimum. The bounds run in three stages, cheapest first: a closed-form bound that disposes of about 99.999% of supports, a capped-simplex QP on the survivors, and an exact enumeration of lot compositions on whatever is left.

The largest cell, 10.27 billion supports, takes 489 s on one GPU. Before we trusted any of it, we checked three things.

• The enumeration is exhaustive: it matches a reference combination generator exactly on four shapes.
• Handed a value above the optimum, it exhibits the better portfolio.
• Running out of budget is reported as unresolved, never as a refutation.

There are two limits.

• No small witness. Unlike the mean-variance certificate, checking this proof means repeating the enumeration.
• It does not scale. At 98 assets and K = 10 there are about 1.5×10¹³ supports, roughly twelve days at the measured 14 million a second.

Mean-CVaR: an impossibility result, then a certificate

Minimise the average of the worst (1 − β) share of scenario losses, less a return reward, holding at most K names. We started with a negative result. The objective is a maximum of functions linear in the weights, and for a linear objective the best K-name portfolio is no better than the best over the convex hull of all K-name portfolios — which, once the caps allow it, is every long-only portfolio. So every bound convex in the weights is capped at the cardinality-free optimum. That covers the LP relaxation, cutting planes, Lagrangian decompositions and Benders cuts alike. On instances with 500 assets and 10,000 scenarios, that ceiling is negative while the best portfolios’ values are positive. Any such bound must therefore report a gap of at least 150% at K = 10 and 276% at K = 20, and the ceiling does not move with K at all.

The route that works gives up on relaxation. With the names fixed, the problem is a linear program, so each support has its own exact dual. That dual gives a lower bound on the support; the minimum over all supports bounds the optimum from below, and a feasible portfolio bounds it from above. Sharing duals across supports failed: 82–97% of supports survived. What worked was noting that the certified bound is a minimum, so effort belongs only on the supports currently holding it down.

These are certified bounds, not proofs to machine precision; the remaining width is the per-support solve’s own convergence. Like index tracking, the method visits every support, so 500 assets remains out of reach, and neither extra problem type yet has an independent checker. The library ships the index-tracking certifier, the CVaR convex-ceiling computation and the per-support CVaR bounds in hc.extras; the enumeration driver behind the table above is still in our research code.

On an idle workstation (a 32-core CPU and an RTX 4080 under WSL2), one factor-tier certificate at n ≈ 1,450 takes 2.6–2.7 s once the solver is compiled and 5.0–5.4 s from a fresh process; at n = 2,901, 4.8 s and 7.8 s. Checking adds 0.9 s in Python and 2.0 s in C at the larger size. Certifying a portfolio is a matter of seconds, warm or cold.

It was not always. Our first full proof at n = 2,901 took 100 s, less than a fifth of it optimisation: 29 s searching for a diagonal split the proposition above says must be zero (a randomised rank test now settles it in 0.1 s), 26 s on an eigenvalue whose only job was to size a safety margin, and about 26 s of just-in-time compilation that can be cached.

Several results that looked like findings about optimisation turned out to be about measurement. Each changed a conclusion.

heurics-cert 1.0.0 is open source under the Apache License 2.0. Loading and checking certificates needs only NumPy; producing them uses JAX. The risk models above are in hc.models (sample_covariance, factor_model, ledoit_wolf, factor_form), and the index-tracking certifier and mean-CVaR bounds are in hc.extras.

The release ships the prover, the three verifiers, the certificate and tree formats, the adversarial test suite, example certificates, and a Dockerfile that runs every check in a clean container. The equity price data cannot be redistributed, so for the 189 mandates it ships a recipe that rebuilds the panel from the public feed, with fingerprints that confirm a rebuild. Every equity model is built without BLAS, so the same returns give byte-identical models — and fingerprints — on any machine.

github.com/Heurics/heurics-cert

Making the copositive split fast at large n. Its ceiling looks flat in size; the obstacle is solving for G, and a second-order method is the natural next step.

Commercial risk models. One forty-factor fundamental model on one date gave tight certificates. Commercial models carry more factors, and more statistical factors loosened the certificate in our sensitivity tests (from a 0.1% median gap with five factors to about 0.9% with twenty). We have not been able to test one.

The other problem types at scale. Index tracking and mean-CVaR are certified by enumeration, which stops at tens of assets for K = 10. For CVaR we have shown that no convex bound can help, so scaling it needs branching with a dual re-optimised at every node. Both also need independent checkers.

The last eight optima. On eight factor-tier mandates the search narrows the gap to at most 0.10% but does not finish.

Independent review of the checkers. All three verifiers have one author. They share no code, but only an outside reader closes the risk of a shared misreading of the mathematics. They are short on purpose and the formats are documented. Independent review is welcome.

Limits. A certificate proves a bound for the model it carries: turnover limits, tax lots or sector bands make a different model. Gap statistics depend on the quality of the best known portfolio. Timings come from one workstation. The exact tier stops near n = 200.