Skip to content

Semidefinite_Programs

Kiguli edited this page Oct 6, 2026 · 3 revisions

Semidefinite Programs

Scen-Opt supports semidefinite programs (SDPs) in inequality form: the scenario constraints and the hard constraints are linear matrix inequalities (LMIs) in a vector of scalar decision variables, and the program has no equality constraints. As with LPs and QPs, the scenario constraints are built from samples $\delta_1, \dots, \delta_N$ of an uncertain parameter $\delta \in \mathbb{R}^q$. After solving, the tool reports the complexity $k$ and the risk bounds $[\underline{\varepsilon}, \overline{\varepsilon}]$; see Scenario Theory.

Formulation

The most general SDP is the relaxed–regularized SDP:

$$\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 & F_0(\delta_i) + \sum_{j=1}^{d} x_j F_j(\delta_i) \preceq \zeta_i\, \mathbb{I}_{\mathfrak{m}}, \qquad i = 1, \dots, N, \\\ & E_0 + \sum_{j=1}^{d} x_j E_j \preceq 0, \end{aligned}$$

with $\tau \ge 0$ and $\rho \ge 0$. Here $M \preceq 0$ means that the symmetric matrix $M$ is negative semidefinite, and $\mathbb{I}_{\mathfrak{m}}$ is the $\mathfrak{m} \times \mathfrak{m}$ identity.

Symbol Meaning
$x = [x_1, \dots, x_d]^\top$ vector of scalar decision variables
$Q \in \mathbb{R}^{d \times d}$ symmetric positive semidefinite matrix of the quadratic cost
$c \in \mathbb{R}^d$ linear cost vector
$F_0(\delta), \dots, F_d(\delta) \in \mathbb{R}^{\mathfrak{m} \times \mathfrak{m}}$ symmetric matrices defining the scenario LMI: one $\mathfrak{m} \times \mathfrak{m}$ LMI per scenario
$E_0, \dots, E_d \in \mathbb{R}^{\mathfrak{n} \times \mathfrak{n}}$ symmetric matrices of the hard LMI, which does not depend on $\delta$ ($\mathfrak{n} = 0$ means none)
$\zeta_i \ge 0$ scalar slack of scenario $i$, which shifts the whole LMI by $\zeta_i \mathbb{I}_{\mathfrak{m}}$
$\rho$; $\tau$, $\bar{x}$, $p$ relaxation penalty; regularization weight, reference point and norm

A matrix variable, such as a Lyapunov matrix $P$, is written in terms of $x$ through a basis: $P = \sum_j x_j \Pi_j$ (see the example below). Each scenario LMI counts as one constraint for the complexity $k$. The hard LMI is never part of the support list.

The four variants

Variant option value Objective Scenario LMIs
Robust robust $c^\top x + \tfrac12 x^\top Q x$ $F_0(\delta_i) + \sum_j x_j F_j(\delta_i) \preceq 0$
Regularized regularization $c^\top x + \tfrac12 x^\top Q x + \tau \lVert x - \bar{x} \rVert_p$ $F_0(\delta_i) + \sum_j x_j F_j(\delta_i) \preceq 0$
Relaxed relaxation $c^\top x + \tfrac12 x^\top Q x + \rho \sum_i \zeta_i$ $F_0(\delta_i) + \sum_j x_j F_j(\delta_i) \preceq \zeta_i \mathbb{I}_{\mathfrak{m}}$, $\zeta_i \ge 0$
Relaxed–regularized regularization-relaxation all four terms $F_0(\delta_i) + \sum_j x_j F_j(\delta_i) \preceq \zeta_i \mathbb{I}_{\mathfrak{m}}$, $\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$.

The hard LMI is present in all four. How to select a variant:

  • Web interface: use the Select Option dropdown. The server ignores the parameters of the variants that are not selected.
  • Python: solve_sdp adds slacks 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 and QPs, the term is the vector $p$-norm of $x - \bar{x}$: any number $p \ge 1$ (for example 1, 1.5, 2 or 3), inf, or fro (the same as $p = 2$). It is added only when $\tau > 0$. Orders other than 1, 2 and inf need a solver with power-cone support, such as MOSEK, CLARABEL or SCS.

Inputs

Symbol Shape Web interface Program JSON key solve_sdp argument
$d$ integer $\ge 1$ dimension field d n_x taken from the length of c
$\mathfrak{m}$ LMI size of the $F_j$ dimension field 𝔪 lmi_size inferred from F_d
$\mathfrak{n}$ LMI size of the $E_j$ (0 = none) dimension field 𝔫 lmi_e_size inferred from E
$c$ $d \times 1$ Edit c c c: (d, 1) or (d,)
$Q$ $d \times d$, no delta Edit Q (optional) Q Q: (d, d), or None/empty for no quadratic term
$F_0(\delta), \dots, F_d(\delta)$ $d+1$ symmetric $\mathfrak{m} \times \mathfrak{m}$ matrices, may use delta[k] Edit F₀,…,d(δ) (symbolic mode) F_d F_d: callable, delta -> {"0": F0, ..., "d": Fd}
$E_0, \dots, E_d$ $d+1$ symmetric $\mathfrak{n} \times \mathfrak{n}$ matrices, no delta Edit E₀,…,d (enabled when 𝔫 > 0) E E: dict {"0": E0, ..., "d": Ed}, or 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 as a 1-D array or a column (default 0)
$p$ see above 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 (d+1)\mathfrak{m}^2$ (numeric) Scenarios upload separate file deltas: 2-D (N, q)

Matrix collections F and E

Both collections are dictionaries keyed by "0", "1", …, "d" (in Python, the integers 0, 1, …, d work too):

  • key "0" holds the constant term ($F_0$ or $E_0$);
  • key "j" holds the matrix that multiplies $x_j$ (the $j$-th entry of $x$, x[j-1] in Python).

In the web interface, the Edit F₀,…,d(δ) and Edit E₀,…,d buttons open a collection editor. It has one button per matrix, plus a box where you can paste the whole collection as JSON. The interface and the program loader require all $d+1$ matrices, each of the size set by 𝔪 (or 𝔫).

  • Symmetry is checked exactly. solve_sdp raises an AssertionError unless every evaluated $F_j(\delta_i)$ and every $E_j$ equals its transpose entry by entry. Write mirrored entries with identical expressions. For example, both off-diagonal entries of $F_3$ in the example below are "0.3*delta[0]-2".
  • E has no delta. Its entries are numbers or constant expressions, so 0.5, "0.5" and "1/2" are all accepted (the matrix editor accepts plain numbers only).
  • F entries. Entries of $F_j(\delta)$ follow the expression rules described for LPs: numbers, delta[0], …, delta[q-1] (indexed from 0), + - * / // % ** and ^, parentheses, functions such as cos(delta[0]) or math.cos(delta[0]), and the constants pi and e; see Linear Programs.
  • JSON entries. In the program JSON, entries of c, Q, F_d and E may be plain numbers or quoted strings; expressions must be quoted. The benchmark files quote every entry.

Requirements on Q

$Q$ must be a symmetric positive semidefinite $d \times d$ matrix. solve_sdp checks symmetry and non-negative eigenvalues with the same small tolerance as solve_qp (see Quadratic Programs), so rank-deficient matrices are accepted.

$Q$ is optional on the SDP tab: an empty $Q$ (or None in Python) means there is no quadratic term, which is the same as a $d \times d$ zero matrix. $d$ is read from the length of $c$.

Notes for Python callers

  • solve_sdp(deltas, F_d, E, 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).
  • Unlike solve_lp/solve_qp, the solution x is a 1-D array of shape (d,).
  • Pass x_ref as a scalar or with d entries, as a 1-D array or a (d, 1) column; any other size raises ValueError.
  • deltas should be 2-D, with shape (N, q). An empty deltas raises ValueError.
  • cost is the value of the full objective, including the $\tau$ and $\rho$ terms. The web interface reports $c^\top x^\ast + \tfrac12 x^{\ast\top} Q x^\ast$.

See SDP.py for implementation details.

Symbolic and numeric input modes

  • Symbolic (default). The $F_j(\delta)$ are entered as expressions. The scenario file is $N \times q$, with one scenario per row.
  • Numeric. The Edit F₀,…,d(δ) button is hidden. Row $i$ of the scenario file holds the evaluated matrices $F_{0,i}, F_{1,i}, \dots, F_{d,i}$, where $F_{j,i} = F_j(\delta_i)$. Each matrix is flattened row by row, and the matrices are concatenated in order $j = 0, \dots, d$. A row therefore has $(d+1) \mathfrak{m}^2$ numbers; the server checks this and reports Each scenario row must have ... elements otherwise. $d$ and $\mathfrak{m}$ must be given (c or n_x, and lmi_size). $E$, $c$ and $Q$ are still entered in the interface or the JSON, because they do not depend on the scenario.

Parameter sweeps

Comma-separated values of $\rho$ and $\tau$ give one solve per $(\tau, \rho)$ pair. Each solve gets its own complexity and risk bounds, computed with the β you entered. The reported overall confidence is $1 - \beta n_\tau n_\rho$ (union bound). See Linear Programs and In-built Graphics.

Program JSON

Key Required Content
type yes "SDP"
mode no "symbolic" (default) or "numeric"
option no "robust", "regularization", "relaxation" or "regularization-relaxation". If absent, it is inferred from the first value of rho and tau, as for LPs.
c yes $d \times 1$ array of numbers or strings
Q no $d \times d$ array of numbers or strings; omit it (or leave it empty) for no quadratic term
F_d symbolic mode object {"0": ..., "d": ...} of $\mathfrak{m} \times \mathfrak{m}$ arrays of numbers or expression strings
E no object {"0": ..., "d": ...} of $\mathfrak{n} \times \mathfrak{n}$ arrays of numbers or constant expressions
n_x, lmi_size, lmi_e_size lmi_size in numeric mode integers $d$, $\mathfrak{m}$, $\mathfrak{n}$. If c, F_d["0"] or E["0"] are present, the dimensions are taken from their shapes.
rho, tau, confidence, p no (but β must be set before solving) strings, e.g. "1.0", "0", "1e-06", "2"
x_ref no $d \times 1$ array

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

Worked example: LPV stability via a common Lyapunov function

Benchmark: benchmarks/SDP_LPV_stability_3d.

Consider the linear parameter-varying (LPV) system $\dot{\xi} = A(\delta) \xi$ with

$$A(\delta) = \underbrace{\begin{bmatrix} 0 & 1 \\ -2 & -1 \end{bmatrix}}_{A_0} + \delta \underbrace{\begin{bmatrix} 0 & 0 \\ 0.3 & 0.1 \end{bmatrix}}_{A_1},$$

where $\delta$ is drawn uniformly from $[-0.22, 1]$; the benchmark uses $N = 100$ samples. We look for a common quadratic Lyapunov function $V(\xi) = \xi^\top P \xi$, with $P \succeq 0$ and $A(\delta)^\top P + P A(\delta) \preceq 0$.

The decision vector holds the three free entries of the symmetric $2 \times 2$ matrix $P$: $x = [p_{11},\ p_{12},\ p_{22}]^\top$. With the basis

$$\Pi_1 = \begin{bmatrix} 1 & 0 \\ 0 & 0 \end{bmatrix}, \quad \Pi_2 = \begin{bmatrix} 0 & 1 \\ 1 & 0 \end{bmatrix}, \quad \Pi_3 = \begin{bmatrix} 0 & 0 \\ 0 & 1 \end{bmatrix}, \qquad P = \sum_{j=1}^{3} x_j \Pi_j,$$

the program data are:

  • Scenario LMI ($\mathfrak{m} = 2$): $F_0(\delta) = 0$ and $F_j(\delta) = A(\delta)^\top \Pi_j + \Pi_j A(\delta)$ for $j = 1, 2, 3$.
  • Hard LMI $P \succeq 0$ ($\mathfrak{n} = 2$): $E_0 = 0$ and $E_j = -\Pi_j$.
  • Objective: $c = [-1,\ 0,\ -1]^\top$, which maximises $\mathrm{trace}(P)$, and $Q = 0.1 \mathbb{I}_3$. The scenario constraints are homogeneous in $P$, so without the quadratic term the trace could grow without bound.
  • Variant: relaxed, with $\rho = 1$, $\tau = 0$ and $\beta = 10^{-6}$.

Program file

data/program_symbolic.json (whitespace reformatted). It has no option key, so "rho": "1.0" makes the loader select Relaxation.

{
  "type": "SDP",
  "c": [["-1"], ["0"], ["-1"]],
  "Q": [["0.1", "0", "0"],
        ["0", "0.1", "0"],
        ["0", "0", "0.1"]],
  "F_d": {
    "0": [["0", "0"], ["0", "0"]],
    "1": [["0", "1"], ["1", "0"]],
    "2": [["-4+0.6*delta[0]", "-1+0.1*delta[0]"], ["-1+0.1*delta[0]", "2"]],
    "3": [["0", "0.3*delta[0]-2"], ["0.3*delta[0]-2", "-2+0.2*delta[0]"]]
  },
  "E": {
    "0": [["0", "0"], ["0", "0"]],
    "1": [["-1", "0"], ["0", "0"]],
    "2": [["0", "-1"], ["-1", "0"]],
    "3": [["0", "0"], ["0", "-1"]]
  },
  "rho": "1.0",
  "tau": "0",
  "confidence": "1e-06",
  "mode": "symbolic",
  "n_x": 3,
  "lmi_size": 2,
  "lmi_e_size": 2
}

data/scenarios.csv has 100 rows with one value of $\delta$ each (delta[0]).

In the web interface

  1. Upload Program with data/program_symbolic.json. This opens the SDP tab with Relaxation, $\rho = 1.0$, $d = 3$, 𝔪 = 2 and 𝔫 = 2.
  2. Upload data/scenarios.csv under Scenarios.
  3. Choose a solver (CLARABEL below) and press Solve.

Results table:

  • optimal cost −9.8439;
  • optimal x = (10.6747, 1.3008, 9.0131);
  • complexity 1;
  • risk bounds [0.0000, 0.1858];
  • confidence 0.999999;
  • no degeneracy;
  • total constraints 101 (the 100 scenario LMIs plus the hard LMI).

In Python

Run from the root of the Scen-Opt repository:

import json
import numpy as np
from src.SDP import solve_sdp
from src.Risk import quantify_risk
from src.parsing import generate_tensor_function, generate_tensor, generate_matrix

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

F_d = generate_tensor_function(json.dumps(prog["F_d"]))  # callable: delta -> {"0": (2,2), ..., "3": (2,2)}
E = generate_tensor(json.dumps(prog["E"]))               # {"0": (2,2), ..., "3": (2,2)}
c = generate_matrix(json.dumps(prog["c"]))               # (3, 1) array
Q = generate_matrix(json.dumps(prog["Q"]))               # (3, 3) array, 0.1 * I

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

# Relaxed SDP: rho = 1, tau = 0
x, zeta, cost, N, k, constraints, degenerate = solve_sdp(
    deltas, F_d, E, c, Q, rho=1.0, solver="CLARABEL")
eps_L, eps_U = quantify_risk(k, N, beta=1e-6)

p11, p12, p22 = x
print("x* = [p11, p12, p22]:", x)
print("objective incl. rho*sum(zeta):", cost)
print("c'x + x'Qx/2:", c.ravel() @ x + 0.5 * x @ Q @ x)
print("largest slack zeta_i:", zeta.max())
print("N =", N, "| k =", k, "| degenerate:", degenerate)
print("risk bounds:", eps_L, eps_U)

Output (CVXPY 1.7.2, CLARABEL):

x* = [p11, p12, p22]: [10.67473652  1.30076817  9.01314939]
objective incl. rho*sum(zeta): -9.843942911792748
c'x + x'Qx/2: -9.843942933747922
largest slack zeta_i: 2.5493322288272623e-10
N = 100 | k = 1 | degenerate: False
risk bounds: 0.0 0.18583116631314622

Reading the result

The Lyapunov matrix found is

$$P^\ast \approx \begin{bmatrix} 10.6747 & 1.3008 \\ 1.3008 & 9.0131 \end{bmatrix},$$

with eigenvalues of about 8.30 and 11.39.

  • MOSEK: the stored MOSEK result in results/metrics.json is $x^\ast = [10.6745,\ 1.3010,\ 9.0133]$.
  • Slacks: all slacks are zero up to solver accuracy. With $\rho = 1$, the relaxed solution therefore satisfies every sampled LMI, and it coincides with the robust solution.
  • Complexity: $k = 1$. Only the LMI of the smallest sample, $\delta \approx -0.187$, is active (its largest eigenvalue is about $-4 \times 10^{-8}$). For every other sample, the largest eigenvalue of $A(\delta_i)^\top P^\ast + P^\ast A(\delta_i)$ is $-0.0245$ or lower.
  • Risk: with confidence at least $1 - 10^{-6}$, the probability that a new $\delta$ from the same distribution violates $A(\delta)^\top P^\ast + P^\ast A(\delta) \preceq 0$ is at most $\overline{\varepsilon} = 0.1858$.

Sweeping $\rho$. Entering ρ = 0.01, 0.1, 1 in the web interface gives:

$\rho$ $x^\ast$ scenarios with $\zeta_i > 0$ $k$ $[\underline{\varepsilon}, \overline{\varepsilon}]$
0.01 (10.6536, 1.2776, 9.0455) 7 8 [0, 0.3137]
0.1 (10.6747, 1.3009, 9.0132) 0 1 [0, 0.1858]
1 (10.6747, 1.3008, 9.0131) 0 1 [0, 0.1858]

The overall confidence is reported as $1 - 3 \times 10^{-6} = 0.999997$.

Numeric mode

For numeric mode, use data/program_numeric.json, which has no F_d and declares "lmi_size": 2, with data/scenarios_numeric.csv. Each row has $(3+1) \cdot 2^2 = 16$ numbers, $[F_{0,i} \mid F_{1,i} \mid F_{2,i} \mid F_{3,i}]$, each matrix flattened row by row. For $\delta_1 \approx -0.18687$, the row is

0, 0, 0, 0,   0, 1, 1, 0,   -4.1121, -1.0187, -1.0187, 2,   0, -2.0561, -2.0561, -2.0374

(rounded here for display). The results are identical to symbolic mode.

Other SDP benchmarks

  • SDP_Quadratic_Stability_6d: common Lyapunov matrix for a 3-state coupled oscillator with two uncertain parameters ($d = 6$, $N = 500$, hard constraint $P \succeq 0.01 \mathbb{I}$).
  • SDP_covariance_wine_15d: minimum-trace covariance matrix that dominates subsample covariances of the UCI Wine data, 5 features ($d = 15$, $N = 178$).
  • SDP_ellipsoid_enclosure_15d: tightest ellipsoid, in the trace sense, containing the UCI Breast Cancer Wisconsin data, 5 features ($d = 15$, $N = 569$, scalar scenario LMIs with $\mathfrak{m} = 1$).

The numeric-mode files of the LPV example are in its own folder (data/program_numeric.json and data/scenarios_numeric.csv, described above).

See also

Using The Tool

Home

Starting Up

Convex Optimization

Graphical Results

Getting into the Weeds

Scenario Theory

Implementation Details

Clone this wiki locally