In a two-good space, initially the consumer maximizes $U(x,z)\;\; s.t. \;\;p_xx+p_zz =I$ and we assume it obtains the solution $(x^*, z^*)$ as a function of prices and income.
In the constrained case, the consumer will either choose $(0, \tilde z)$ or $(x^*+\epsilon, z'$), for some $\epsilon >0 $ always exhausting its budget, so in particular, $\tilde z = I/p_z$.
In order for the consumer to still choose to buy a strictly positive quantity of $x$, it must be the case that
$$U(x^*+\epsilon, z') > U(0, \tilde z)$$
Apply a first order approximation around $(x^*, z^*)$ without ignoring the remainders, we want
\begin{align}
U(x^*, z^*) + U_x(x^*)\cdot \epsilon + U_z(z^*)(z'-z^*) + R_{\epsilon} &\\> U(x^*, z^*) + U_x(x^*)(-x^*) + U_z(z^*)(\tilde z-z^*) + R_0
\end{align}
Simplify and re-arrange, we want
$$U_x(x^*)(x^*+\epsilon) + R_{\epsilon} > U_z(z^*)(\tilde z-z') + R_0 $$
We know that from the unconstrained optimization, $U_x(x^*)/U_z(z^*) = p_x/p_z$ so
$$\frac {p_x}{p_z}(x^*+\epsilon) + \frac {R_{\epsilon}}{U_z(z^*)} > \left(\frac{I}{p_z}-z'\right) + \frac {R_0}{U_z(z^*)} $$
Multiply throughout by $p_z$,
$$p_x(x^*+\epsilon) + p_z\frac {R_{\epsilon}}{U_z(z^*)} > I - p_zz' + p_z\frac {R_0}{U_z(z^*)} $$
but $p_x(x^*+\epsilon) + p_zz' = I \implies p_x(x^*+\epsilon) = I-p_zz'$ so we are left with the requirement that (ignoring positive terms)
$$R_{\epsilon} > R_0 $$
in order for the consumer to choose $x^*+ \epsilon$ and not $0$ for $x$.
Note that the above take into account also the signs of the remainders, it is not just about their absolute magnitudes.
Now let's go back to our first-order expansions. We know that both candidate bundles yield utilities lower than $U(x^*, z^*)$, because they were feasible in the unconstrained case, and they weren't chosen.
Looking at the expansion of $U(0, \tilde z)$ we then conclude that we have
$$U_x(x^*)(-x^*) + U_z(z^*)(\tilde z-z^*) + R_0 < 0 $$
$$\implies U_z(z^*)\cdot \Big[(U_x(x^*)/U_z(z^*))\cdot(-x^*) + \tilde z-z^*\Big] + R_0 < 0$$
$$\implies \frac {U_z(z^*)}{p_z}\cdot \Big[-p_xx^* + p_z\tilde z-p_zz^*\Big] + R_0 < 0$$
But $-p_xx^* -p_zz^* = -I$ and $p_z\tilde z =I$ so the term in brackets is zero. Therefore we conclude that
$$R_0 <0$$
Looking now at the expansion of $U(x^*+\epsilon, z')$, we know we have
$$U_x(x^*)\cdot \epsilon + U_z(z^*)(z'-z^*) + R_{\epsilon} < 0$$
Performing the same manipuations as before we obtain here too that
$$R_{\epsilon} < 0$$
So the condition to buy $x^*+\epsilon$ can be re-written as
$$|R_{\epsilon}| < |R_0|$$
This formalizes somewhat the notion that if $\epsilon$ is "sufficiently small", $R_{\epsilon}$ will be smaller in absolute terms than $R_0$, since the approximation to the same function will be "better", and so we will observe $x^*+\epsilon$ and not $0$. But it also tells us what the graphs in the other answer told us too, that there is not a single general answer to the matter.
