How can I find the largest axis-parallel square that can be inscribed into a linearly transformed unit square? The initial square is centered about the origin and each vertex $\boldsymbol v_i \in \{-1,1\}^2 \subset \mathbb R^2$ is transformed by $\boldsymbol M\boldsymbol v_i$. The transformation matrix $\boldsymbol M$ is an arbitrary $2\times 2$ matrix. Now, I want to find the largest axis-parallel and origin-centered square that can be inscribed into the transformed square, e.g. expressed by a scaling transform of the initial square.
My first approach was to state this as a maximization problem to find the scaling factor $\lambda^{\ast}$ like the following:
$$\lambda^{\ast} = \arg\max_{\lambda} \lambda$$
$$s.t. \quad \boldsymbol M(\operatorname{rot_{-90^{\circ}}}(\boldsymbol v_{i+1 \operatorname{mod} 4} - \boldsymbol v_{i}))\cdot \boldsymbol (\lambda\boldsymbol v_k - \boldsymbol M\boldsymbol v_i) \leqslant 0 \quad \text{for $i,k \in \{0,1,2,3\}$},$$
where the indices enumerate the vertices in counter-clockwise order. But I am not sure how to continue from here.
 A: This is a problem whose answer is geometrically beautiful but hard to formalise.

If $\mathbf M=\begin{bmatrix}a&b\\c&d\end{bmatrix}$, three of the points in the transformed square will be $A=(a+b,c+d)$, $B=(-a+b,-c+d)$ and $C=(a-b,c-d)$. The line through $A$ and $B$ has equation $l_1:ay-cx=a(-c+d)-c(-a+b)$; that through $A$ and $C$, $l_2:by-dx=b(c-d)-d(a-b)$.
Now $l_1$ and $l_2$ have to be intersected with the lines $d_+:x+y=0$ and $d_-:x-y=0$ (i.e. $y=\pm x$), making for four intersection points. Any intersections at infinity (parallel lines) should be disregarded; the finite intersection points have the form $(s_1a,s_2a)$ where $a\ge0$ and $s_1$ and $s_2$ are signs. The smallest $a$ among the intersections, $a^*$, is half the side length of the largest axis-parallel, origin-centred square $S^*$ that fits in the transformed square – the coordinates are thus $\{\pm a^*,\pm a^*\}$.
Why is this so? The vertices of $S^*$ lie on $y=\pm x$ as well, and both it and the transformed square are convex. Therefore, if rays are cast out along $d_\pm$ from the origin, the shortest distance among them will be half the diagonal of $S^*$. The above algorithm is a formalisation of this intuition.
