Notation. $H^k(\mathbb R^n)$ denotes the Sobolev space $W^{k,2}(\mathbb R^n)$. The Sobolev space $W^{k,p}(\mathbb R^n)$ is the Banach space of functions in $L^p$ with all derivatives up to order $k$ in $L^p$.
My question is about non-uniqueness of weak solutions of an elliptic PDE.
The PDE
Consider a very simple semilinear elliptic PDE: $$ \left\{\begin{aligned} &-\Delta u=\partial_{x_1}(u^2+f)& &\text{in } \mathbb R^n,\\ & \,\,\,u(x)\to 0\,,& &x\to\infty. \end{aligned}\right. $$ where $x=(x_1,...,x_n)$. The boundary conditions are not rigorous, it is just to describe the kind of solutions we want. The critical Lebesgue spaces with respect to the natural scaling symmetry are: $$ u\in L^{n}(\mathbb R^n), f\in L^{n/2}(\mathbb R^n). $$
Strong solutions
Assume $n=3$ to fix ideas. By using the fixed point theorem, it is possible to show that if $f\in L^{3/2}(\mathbb R^3)$ is small enough, then there exists a unique small solution $u\in L^3(\mathbb R^3)$, which moreover satisfies $$ \|\nabla u\|_{L^{3/2}(\mathbb R^3)}\lesssim \|f\|_{L^{3/2}(\mathbb R^3)} $$ (note that $\dot W^{1,3/2}(\mathbb R^3)\hookrightarrow L^3(\mathbb R^3)$).
Weak solutions
If I instead consider a scaling-subcritical datum, like $f\in L^2(\mathbb R^3)$, I cannot use the fixed point theorem. But notice that multiplying the equation by $u$ and integrating by parts, we obtain the a priori bound
$$ \|\nabla u\|_{L^2(\mathbb R^3)}\lesssim \|f\|_{L^2(\mathbb R^3)}. $$
Using this a-priori bound (please correct me if I am wrong*) it is possible to prove the existence of weak solutions $u\in \dot H^1(\mathbb R^3)$ for any $f\in L^2(\mathbb R^3)$ (we have $\dot H^1(\mathbb R^3)\hookrightarrow L^6(\mathbb R^3)$).
In this case, estimates for the difference of two solutions do not seem to exist, so the uniqueness cannot be proved by standard arguments.
My question
My problem is the following. Is the weak solution unique if $f\in L^2(\mathbb R^3)$? If no, how can one find a counterexample?
In a case where we have a priori bounds like in the PDE I wrote, do we have non-uniqueness in case of subcritical data? Note that this is not a regularity prolem: data in $L^2$ are "more regular" than data in $L^{3/2}$. The problem is that data grow more slowly at infinity than data in the critical space. I could ask the same question for data $f\in H^{100}(\mathbb R^3)$, and the answer would still be non trivial. Also, I am reasonably sure that a unique continuation property exists for this kind of equation, so this problem is extremely different from those of non-uniqueness of weak solutions of evolution PDE at low regularity (like, e.g., the famous results for Euler equation).
I feel it is likely that there are data that lead to two different solutions, but I know no technique that would do that or examples like that in the literature. If you know of articles about anything that resembles this problem, any reference is welcome.
*I think there should be more than one way of doing that. The way I thought would be to approximate the Cauchy problem with that of the same PDE on a half plane (with Dirichlet boundary conditions), where the "edge" of the plane goes to infinity. On a half plane, you can look at the PDE as an evolution equation via the Poisson kernel, so you can show that a solution exists, and then you can use the uniform bound and weak compactness to show the existence of the solution of the original problem... Probably not the easiest way. I have to say I am not an elliptic person, I am mostly dispersive, sometimes parabolic... In any case, feel free to comment on this point.