Underconstrained linear system, finding solution closest to $w_0\ne 0$?

Question

How do I find solution of an under-constrained linear system in $\mathbb{R}^d$ that is closest to given reference $w_0$ in Euclidean distance?

Below is a solution that uses QuadraticOptimization, however that requires writing out the formulas symbolically. For $d=1000$, there are 1 million symbolic variables, is there a more efficient approach?

norm2[mat_] := Total@Flatten[mat*mat];
modifiedLeastSquares[m_, b_, w0_] := Module[{xx, x, eqs, d},
   d = Last@Dimensions[b];
   x = Array[xx, {d, d}];
   eqs = Flatten@MapThread[#1 == #2 &, {m . x, b}, 2];
   x /. QuadraticOptimization[norm2[x - w0], eqs, Flatten@x]
   ];

SeedRandom[1];
bs = 4;
d = 8;
m = RandomReal[{-1, 1}, {bs, d}];
b = RandomReal[{-1, 1}, {bs, d}];
w0 = ConstantArray[0, {d, d}];
mySol = modifiedLeastSquares[m, b, w0];
refSol = LeastSquares[m, b];
Print["Difference from built-in LeastSquares: ", 
 mySol - refSol // Norm]

Answer 1

Your problem

• Find $x$ such that $m.x=b$ and $||x-w_{0}||^2\to $ min

is equivalent to (by change of variables $y=x+w_{0}$)

• Find $y$ such that $m.y=b-m.w_{0}$ and $||y||^2\to $ min

So, the simplest answer is:

LeastSquares[m, b - m . w0] + w0

This solves the case of $d=1000$ and $bs=500$ within 1.2 seconds on my PC.