Minimizing a squared error loss with nuclear norm constraint I have to minimize an objective function of the the form :

$||X_{s} - Y_{s}D_{s}||_{F}^{2} + (||D_{s}||_{*}^{2} - 1)$

where $||\cdot||_{F}$ denotes the Frobenius norm and $||\cdot||_{*}$ denotes the nuclear norm.
In the above equation, $X_{s}$ and $Y_{s}$ are data matrices and an optimal value of $D_{s}$ needs to be found out so that it not only has a low value for the nuclear norm but also is able to reconstruct $X_{s}$ given $Y_{s}$.
From the literature available for the nuclear norm, I don't think there is a closed form solution for the above optimization problem ? Also, I am not trying interior point methods since, in the literature, they are reported to have given sub-optimal results.
I was trying to optimize the above function using CVXPY library in Python but unfortunately the code used to collapse  in the middle citing segmentation error.
A similar question was posted by me on How to minimize the sum of Frobenius norm and Nuclear norm
Could you please help me with some alternate way (preferably using any relevant Python based library) to get the optimal solution for the above problem ?
 A: The proximal forward-backward splitting algorithm of Combettes and Wajs might be what you want here.  PyUNLocBoX: Optimization by Proximal Splitting is a  Python package for it (Disclaimer: I have only ever written my own code for this algorithm, so I have never used this package), and documentation on the forward-backward algorithm for the package is here. You will probably have to adapt their function slightly, as it is written for vector variables, not matrix variables.
The proximal forward-backward splitting algorithm is an iterative algorithm for finding the minimum of a function $h(x) = f(x) + g(x)$, where both $f$ and $g$ are proper convex functions that are lower-semicontinuous, and where $f$ has an $L$-Lipschitz continuous gradient. The function $g$ does not have to be differentiable, however. So you could write $f(D) = \| X - YD \|_F^2$ and $g(D) = \|D\|_*^2$ (In your objective function, the $-1$ is unnecessary, as it does not affect the argmin. That's one reason why I think you might have meant to say $(\|D\|_* - 1)^2$.
The updates are of the form
$$
    x_{k + 1} = [\text{prox}_{\eta g} \circ (\text{id} - \gamma \nabla f)](x_k)
$$
where $\text{prox}_{\eta g}$ is the proximal operator associated with the function $\eta g$ for some parameter $\eta > 0$, and where $0 < \gamma < \frac{2}{L}$. In effect, the iterative step is really two different iterative steps, one which takes a gradient descent step to try to reduce the differentiable piece $f(x)$, and then a proximal-point (kind of like gradient descent for non-differentiable functions) step which to try to reduce the non-differentiable piece $g(x)$. Here, $\gamma$ is the step size for the first part of the update step, and $\eta$ is the step size for the second part of the update step.
The proximal operator associated with a function $g \colon \mathbb{R}^n \to \mathbb{R}$ is defined as
$$
\text{prox}_{\eta g} (x) = \text{argmin}_{w \in \mathbb{R}^n} \left( g(w) + \frac{1}{2 \eta} \| w - x \|_2^2 \right)
$$
and you can compute this operator for the nuclear norm, the proximal operator associated to the function $\eta \| \cdot \|_* \colon \mathbb{R}^{m \times n} \to \mathbb{R}$ is given by the soft-thresholding operator with parameter $\eta$ applied to the singular values of the argument, i.e.,
$$
\text{prox}_{\eta \| \cdot \|_*} (D) = U \mathcal{S}_\eta (\sigma (D)) V^{\text{T}}
$$
where $\mathcal{S}_\eta$ is the soft-thresholding operator with parameter $\eta$ and svd$ (D) = U \sigma(D) V^{\text{T}}$ is the singular value decomposition of $W$. (This is the part you might have to adapt.)
