Motivation

Many problems in science and engineering can be formulated as a constrained least-squares problem

$$ \begin{aligned} \mathrm{minimize} & \quad | \mathbf{Ax} - \mathbf{b} |^2 \\ \mathrm{subject~to} & \quad \mathbf{Cx} = \mathbf{d}. \end{aligned} $$

This is the case for instance for simple linear optimal control without any bound constraints on the actuation or output signals.

stdlib already provides the unconstrained lstsq solver based on gelsd. LAPACK however also exposes the gglse routines solving precisely the constrained version.

Prior Art

As far as I know, neither Python nor Julia provide an equality constrained least-squares solver out of the box. There are certainly some dedicated packages for that task but the implementations ain't part of numpy/scipy (Python) or LinearAlgebra (Julia).

Additional Information

The calling syntax could be something along the lines of:

x = constrained_lstsq(A, b, C, d[, overwrite_matrices, err])

In contrast to gelsd for unconstrained least-squares, the gglse for equality-constrained least-squares does not provide estimates of the condition number nor of the rank of the matrices. A subroutine version might also be exposed for advanced usage.

Ping: @perazz, @jalvesz, @jvdp1