With the spin-orbit force, we can modify our Woods-Saxon potential to
$$
\hat{u}_{\mathrm{ext}}(r)=-\frac{V_0}{1+\exp{(r-R)/a}}+V_{so}(r)\boldsymbol{ls},
$$
with
$$
V_{so}(r) = V_{so}\frac{1}{r}\frac{d f_{so}(r)}{dr},
$$
where we have
$$
f_{so}(r) = \frac{1}{1+\exp{(r-R_{so})/a_{so}}}.
$$