Quadratic_Programs

Quadratic Programs

A scenario quadratic program (QP) is a scenario LP with a convex quadratic term added to the cost. The constraints are the same linear constraints as in the LP: scenario constraints built from the samples $\delta_1, \dots, \delta_N$ of an uncertain parameter $\delta \in \mathbb{R}^q$, plus optional hard constraints. After solving, Scen-Opt reports the complexity $k$ and the risk bounds $[\underline{\varepsilon}, \overline{\varepsilon}]$; see Scenario Theory.

QPs handle inputs, input modes and parameter sweeps exactly as LPs do. This page restates the essentials and adds what is specific to $Q$.

Formulation

The most general QP is the relaxed–regularized QP:

$$\begin{aligned} \min_{x \in \mathbb{R}^d,\ \zeta_1,\dots,\zeta_N \ge 0} \quad & c^\top x + \tfrac{1}{2} x^\top Q x + \tau \lVert x - \bar{x} \rVert_p + \rho \sum_{i=1}^{N} \zeta_i \\\ \text{subject to} \quad & A(\delta_i)\, x + b(\delta_i) \le \zeta_i, \qquad i = 1, \dots, N, \\\ & G\, x + h \le 0, \end{aligned}$$

with $\tau \ge 0$ and $\rho \ge 0$.

Symbol	Meaning
$x \in \mathbb{R}^d$	decision vector
$Q \in \mathbb{R}^{d \times d}$	symmetric positive semidefinite matrix of the quadratic cost
$c \in \mathbb{R}^d$	linear cost vector
$A(\delta) \in \mathbb{R}^{\mathfrak{m} \times d}$, $b(\delta) \in \mathbb{R}^{\mathfrak{m}}$	scenario constraints ($\mathfrak{m}$ rows per scenario)
$G \in \mathbb{R}^{\mathfrak{n} \times d}$, $h \in \mathbb{R}^{\mathfrak{n}}$	hard constraints ($\mathfrak{n} = 0$ means none); never part of the support list
$\zeta_i \ge 0$	slack of scenario $i$: one scalar, shared by the $\mathfrak{m}$ rows of scenario $i$
$\rho$; $\tau$, $\bar{x}$, $p$	relaxation penalty; regularization weight, reference point and norm

Each scenario counts as one constraint, so the complexity $k$ counts scenarios, not rows.

The four variants

Variant	`option` value	Objective	Scenario constraints
Robust	`robust`	$c^\top x + \tfrac12 x^\top Q x$	$A(\delta_i)x + b(\delta_i) \le 0$
Regularized	`regularization`	$c^\top x + \tfrac12 x^\top Q x + \tau \lVert x - \bar{x} \rVert_p$	$A(\delta_i)x + b(\delta_i) \le 0$
Relaxed	`relaxation`	$c^\top x + \tfrac12 x^\top Q x + \rho \sum_i \zeta_i$	$A(\delta_i)x + b(\delta_i) \le \zeta_i$, $\zeta_i \ge 0$
Relaxed–regularized	`regularization-relaxation`	all four terms	$A(\delta_i)x + b(\delta_i) \le \zeta_i$, $\zeta_i \ge 0$

The variants relate to the general formulation as follows:

Robust: $\tau = 0$, $\rho = 0$ and $\zeta_i = 0$.
Relaxed: $\tau = 0$.
Regularized: $\rho = 0$ and $\zeta_i = 0$.

How to select a variant:

Web interface: use the Select Option dropdown. It shows only the fields the selected variant uses, and the server ignores the parameters of the other variants.
Relaxed options with $\rho = 0$: the web interface keeps the slacks, so violations cost nothing; use $\rho > 0$.
Python: solve_qp adds slack variables when include_slack=True. With the default include_slack=None, it adds them whenever rho != 0.

Norm p. Set it with the web field p, the JSON key p, or norm_type in Python; the default is 2. As for LPs, the term is the vector $p$-norm of $x - \bar{x}$: any $p \ge 1$, inf, or fro (the same as $p = 2$). It is added only when $\tau > 0$. With $p = 1$ or inf a QP solver such as OSQP can be used; other orders need a conic solver (MOSEK, CLARABEL, SCS).

Inputs

Symbol	Shape	Web interface	Program JSON key	`solve_qp` argument
$d$, $\mathfrak{m}$, $\mathfrak{n}$	integers	dimension fields d, 𝔪, 𝔫	`n_x`, `rows_A`, `rows_G`	inferred
$c$	$d \times 1$	Edit c	`c`	`c`: `(d, 1)` or `(d,)`
$Q$	$d \times d$, no `delta`	Edit Q (required)	`Q`	`Q`: `(d, d)`
$A(\delta)$	$\mathfrak{m} \times d$, may use `delta[k]`	Edit A(δ) (symbolic mode)	`A_d`	callable `delta -> (m, d)`
$b(\delta)$	$\mathfrak{m} \times 1$, may use `delta[k]`	Edit b(δ) (symbolic mode)	`b_d`	callable `delta -> (m, 1)` or `(m,)`
$G$, $h$	$\mathfrak{n} \times d$, $\mathfrak{n} \times 1$, numeric	Edit G, Edit h (when 𝔫 > 0)	`G`, `h`	arrays, or `np.array([])` for none
$\tau$, $\rho$	scalars ≥ 0; comma-separated lists in the web/JSON	τ ≥ 0, ρ ≥ 0	`tau`, `rho`	`tau`, `rho` (floats, default 0)
$\bar{x}$	$d \times 1$	Edit x̄	`x_ref`	`x_ref`: scalar or `d` entries (default 0)
$p$	any number $\ge 1$, `inf`, `fro`	p	`p`	`norm_type` (default 2)
$\beta$	$0 < \beta < 1$	Confidence ... where β	`confidence`	`beta` in `quantify_risk(k, N, beta)`
scenarios	$N \times q$ (symbolic) or $N \times (\mathfrak{m}d + \mathfrak{m})$ (numeric)	Scenarios upload	separate file	`deltas`: 2-D `(N, q)`

Requirements on Q

solve_qp checks $Q$ before solving and raises an AssertionError (Q must be symmetric or Q must be positive semi-definite (eigenvalues: ...)) if either check fails:

$Q$ must equal $Q^\top$ up to a small tolerance.
The smallest eigenvalue of $Q$ must be $\ge 0$, again up to a small tolerance. Rank-deficient matrices assembled in floating point, such as $D^\top D$, whose zero eigenvalues come out as about $-10^{-16}$, are accepted.

Both tolerances are $10^{-9}$ times the largest absolute entry of $Q$ (or $10^{-9}$ if all entries are below 1). See QP.py for details.

To drop the quadratic term, use a zero matrix. The QP tab always requires a $d \times d$ matrix $Q$.

Other rules (same as for LPs)

Expressions. Entries of $A(\delta)$ and $b(\delta)$ are numbers or expressions in delta[0], …, delta[q-1], indexed from 0. They may use + - * / // % ** (and ^ for power), parentheses, functions such as sin, exp, sqrt and abs, the constants pi and e, and math. functions such as math.sin(delta[1]); the full list is in Linear Programs. $c$, $Q$, $G$ and $h$ cannot contain delta.
JSON entries. In the program JSON, entries may be plain numbers or quoted strings (0, "1", "-delta[0]*delta[2]"); expressions must be quoted.
Input modes.
- Symbolic: the scenario file is $N \times q$, with one scenario per row.
- Numeric: row $i$ is $A(\delta_i)$ flattened row by row, followed by $b(\delta_i)$, so each row has $\mathfrak{m}d + \mathfrak{m}$ numbers. $d$ and $\mathfrak{m}$ must be given (c or n_x, and rows_A).
Details and an example are in Linear Programs.
Sweeps. Comma-separated values of $\rho$ and $\tau$ give $n_\tau n_\rho$ solves. The reported overall confidence is $1 - \beta n_\tau n_\rho$ (union bound).
Python shapes. deltas should be 2-D, (N, q). b_d may return (m, 1) or (m,) and h may be (n, 1) or (n,); both are reshaped to columns.

solve_qp(deltas, A_d, b_d, G, h, c, Q, tau=0.0, x_ref=np.array([0.0]), rho=0.0, norm_type=2, solver=None, include_slack=None) returns (x, zeta, cost, N, complexity, constraints, degeneracy):

x has shape (d, 1).
cost is the value of the full objective, including the $\tau$ and $\rho$ terms. The web interface instead reports $c^\top x^\ast + \tfrac12 x^{\ast\top} Q x^\ast$.

See QP.py for implementation details.

Program JSON

The program JSON uses the same keys as for LPs (see Linear Programs), with "type": "QP" and the extra required key Q ($d \times d$ array of numbers or strings). The keys are:

type;
mode;
option, which is inferred from rho/tau if absent;
c, Q;
A_d, b_d, required in symbolic mode;
G, h, optional;
n_x, rows_A (required in numeric mode), rows_G;
rho, tau, confidence, p, x_ref.

Optional fields that you do not use can be omitted or left empty ([], "" or null).

Worked example: support vector machine on Iris petal features

Benchmark: benchmarks/QP_Iris_minimal_3d.

The goal is to separate Iris setosa from the other two species using two features, petal length and petal width, with a hard-margin linear support vector machine (SVM). Each of the $N = 150$ scenarios is a triple $\delta_i = (\delta_i[0], \delta_i[1], \delta_i[2])$:

$(\delta_i[0], \delta_i[1])$ are petal length and width in cm;
$\delta_i[2] \in \lbrace -1, +1 \rbrace$ is the label: $-1$ for the 50 setosa samples, $+1$ for the other 100.

Write $z_i = (\delta_i[0], \delta_i[1])$ for the features and $y_i = \delta_i[2]$ for the label. The separating line is $\mathbf{w}^\top z + b = 0$, and the decision vector is $x = [\mathbf{w}_1,\ \mathbf{w}_2,\ b]^\top$. The hard-margin SVM

$$\min_{\mathbf{w}, b}\ \tfrac12 \lVert \mathbf{w} \rVert^2 \quad \text{subject to} \quad y_i\,(\mathbf{w}^\top z_i + b) \ge 1$$

is a robust QP with

$$c = \begin{bmatrix} 0 \\ 0 \\ 0 \end{bmatrix}, \quad Q = \begin{bmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 0 \end{bmatrix}, \quad A(\delta_i) = \begin{bmatrix} -\delta_i[0]\,\delta_i[2] & -\delta_i[1]\,\delta_i[2] & -\delta_i[2] \end{bmatrix}, \quad b(\delta_i) = 1,$$

so each scenario contributes one row ($\mathfrak{m} = 1$). There are no hard constraints, $\tau = \rho = 0$ and $\beta = 10^{-6}$.

Program file

data/program_symbolic.json (whitespace reformatted):

{
  "type": "QP",
  "A_d": [["-delta[0]*delta[2]", "-delta[1]*delta[2]", "-delta[2]"]],
  "b_d": [["1"]],
  "c": [["0"], ["0"], ["0"]],
  "Q": [["1", "0", "0"],
        ["0", "1", "0"],
        ["0", "0", "0"]],
  "rho": "0",
  "tau": "0",
  "confidence": "1e-06",
  "mode": "symbolic",
  "n_x": 3,
  "rows_A": 1,
  "rows_G": 0
}

data/scenarios.csv has 150 rows of the form petal_length, petal_width, label. The first row is 1.4,0.2,-1.

In the web interface

Upload Program with data/program_symbolic.json. This opens the QP tab with Robust, $d = 3$, 𝔪 = 1 and 𝔫 = 0.
Upload data/scenarios.csv under Scenarios.
Choose a solver (CLARABEL below) and press Solve.

Results table:

optimal cost 1.1765;
optimal x = (1.2941, 0.8235, −3.7882);
complexity 2;
risk bounds [0.0000, 0.1441];
confidence 0.999999;
no degeneracy;
total constraints 150.

In Python

Run from the root of the Scen-Opt repository:

import json
import numpy as np
from src.QP import solve_qp
from src.Risk import quantify_risk
from src.parsing import generate_matrix_function, generate_matrix

data = "benchmarks/QP_Iris_minimal_3d/data/"
prog = json.load(open(data + "program_symbolic.json"))

A_d = generate_matrix_function(json.dumps(prog["A_d"]))  # callable: delta -> (1, 3) array
b_d = generate_matrix_function(json.dumps(prog["b_d"]))  # callable: delta -> (1, 1) array
c = generate_matrix(json.dumps(prog["c"]))               # (3, 1) array
Q = generate_matrix(json.dumps(prog["Q"]))               # (3, 3) array, diag(1, 1, 0)
G, h = np.array([]), np.array([])                        # no hard constraints

deltas = np.loadtxt(data + "scenarios.csv", delimiter=",", ndmin=2)  # shape (150, 3)

# Robust QP: tau = 0 and rho = 0 (the defaults)
x, zeta, cost, N, k, constraints, degenerate = solve_qp(
    deltas, A_d, b_d, G, h, c, Q, solver="CLARABEL")
eps_L, eps_U = quantify_risk(k, N, beta=1e-6)

w1, w2, b = x.ravel()
print("x* = [w1, w2, b]:", x.ravel())
print("cost (1/2)||w||^2:", cost, "| margin 2/||w||:", 2 / np.hypot(w1, w2))
print("N =", N, "| k =", k, "| degenerate:", degenerate)
print("risk bounds:", eps_L, eps_U)

Output (CVXPY 1.7.2, CLARABEL):

x* = [w1, w2, b]: [ 1.29411765  0.82352941 -3.7882353 ]
cost (1/2)||w||^2: 1.1764705890324765 | margin 2/||w||: 1.3038404805987853
N = 150 | k = 2 | degenerate: False
risk bounds: 0.0 0.14409952527843412

Reading the result

Classifier: the separating line is $1.2941 z_1 + 0.8235 z_2 = 3.7882$. Points with $\mathbf{w}^\top z + b > 0$ are classified as not setosa. The margin width $2/\lVert \mathbf{w} \rVert$ is 1.3038.
Complexity: two support vectors determine the solution, so $k = 2$.
Risk: with confidence at least $1 - 10^{-6}$, a new iris drawn from the same population violates $y (\mathbf{w}^\top z + b) \ge 1$ (it lies inside the margin or is misclassified) with probability at most $\overline{\varepsilon} = 0.1441$.
MOSEK: the stored MOSEK result in results/metrics.json agrees to three decimals.

If the classes were not linearly separable, the robust QP would be infeasible. The Relaxation option (with $\rho > 0$) then gives the soft-margin SVM, where $\zeta_i$ measures how far sample $i$ violates its margin constraint.

Numeric mode

For numeric mode, use data/program_numeric.json, which has no A_d/b_d and declares "rows_A": 1, with data/scenarios_numeric.csv. Each row is $[A(\delta_i) \mid b(\delta_i)]$, i.e. 4 numbers. The first row is 1.4, 0.2, 1, 1. The results are identical to symbolic mode.

Other QP benchmarks

QP_cbf_robot_navigation_118d: trajectory of a point robot with double-integrator dynamics (20 steps) around a wall whose position is uncertain ($d = 118$, 320 hard constraints, $N = 500$).
QP_radiation_therapy_50d: prostate brachytherapy plan under catheter-placement uncertainty, built on TROTS dose-influence data ($d = 50$, 80 dose constraints per scenario). It is a regularized QP with $\tau = 0.1$. Degeneracy is detected: the published MOSEK result has $k = 5$, and the web interface shows such a result's risk bounds as [N/A, 0.1414].

Navigation Menu

Repository navigation

Quadratic_Programs

Quadratic Programs

Formulation

The four variants

Inputs

Requirements on Q

Other rules (same as for LPs)

Program JSON

Worked example: support vector machine on Iris petal features

Program file

In the web interface

In Python

Reading the result

Numeric mode

Other QP benchmarks

See also

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Using The Tool

Starting Up

Convex Optimization

Graphical Results

Getting into the Weeds

Scenario Theory

Implementation Details

Clone this wiki locally