Harnack's inequality concerns about the non-negative function which is harmonic in $U$, by the choice of $r$, $x$, and $y$, we know that $B(y,r)\subset B(x,2r) \subset U$, hence the inequality holds.
First note that if $\alpha<\frac{n-p}{p}$, then the function $v(x)=\frac{1}{|x|^\alpha}$ is such that $v\in W^{1,p}(U)$. You can prove this by using spherical coordinates. Define $$u_k(x)=\left(\sum _{i=1}^k\frac{1}{2^i|x-r_i|^\alpha}\right)^p $$
Note that:
I - $u_k(x)\leq u_{k+1}(x)$ almost everywhere in $U$.
II - $\sup_{k} \int_U u_k<\infty$
I is immediate. To prove II, note that $$\int_U u_k(x)dx=\int_U\sum_{i=1}^k\left(\frac{1}{2^i|x-r_i|^\alpha}\right)^p\tag{1}$$
We can apply Minkowski Integral Inequality in $(1)$ to get $$\int_U u_k(x)dx\leq \left(\sum_{i=1}^k\left(\int_U\frac{1}{2^{ip}|x-r_i|^{\alpha p}}\right)^{1/p}\right)^p$$
Because $v\in L^p(U)$, we conclude from $(1)$ that $$\sup_{k}\int_U u_k(x)dx<\infty$$
From I and II, we can apply the monotone convergence theorem to conclude that $u$ is well defined and $u_k\to u^p$ in $L^1(U)$, which is the same to say that $u_k^{1/p}\to u$ in $L^p(U)$.
Now define $$g(x)=\sum_{i=1}^\infty \nabla \left(\frac{1}{2^i|x-r_i|^\alpha}\right)$$
By using a similar argument as above, you can prove that the sequence of partial sums $(g_k)$ of $g$, converge to $g$ in $L^p(U)$, which implies that $g\in L^p (U)$. Now note that $$\int_U u_k(x)\frac{\partial\varphi(x)}{\partial x_i}dx=-\int_U \frac{\partial g_{i,k}}{\partial x_i}\varphi(x),\ \forall\ \varphi\in C_0^{\infty}(U)\tag{2}$$
where $g_{i,k}=\frac{\partial g_k}{\partial x_i}$. By taking limit in both sides of $(2)$, we conclude that
$$\int_U u(x)\frac{\partial\varphi(x)}{\partial x_i}dx=-\int_U \frac{\partial g_i}{\partial x_i}\varphi(x),\ \forall\ \varphi\in C_0^{\infty}(U)\tag{3}$$
where $g_i$ is the $i$-th entry of $g$. From $(3)$ we conclude that $u\in W^{1,p}(U)$ and $\nabla u=g$.
To answer your third question, note that $u(r_k)=\infty$ for every $k$, and because $(r_k)$ is dense in $U$, we conclude that $u$ is unbounded in each open set contained in $U$.
Best Answer
The first equation is a special case of the second one. Pick any x in $B⁰(0, r)$. Clearly $u(x) \leq \sup_V u(\cdot)$ holds and, due to Harnack’s inequality, $u(x) \leq \sup_V u(\cdot) \leq C \inf_V u(\cdot)$ is true. Since $\inf_V u(\cdot) \leq u(y)$ holds for any $y$, we can conclude $u(x) \leq C u(0)$. If we interchange the roles of $x$ and $0$, we can $u(0) \leq C u(x)$ or $\frac{1}{C} u(0) \leq u(x)$. Putting this together, we get $\frac{1}{C} u(0) \leq u(x) \leq C u(0)$. This is the first statement.
Addition: This is not exactly the first statement since the constants aren’t reciprocal to each other. The first equation can therefore be regarded as a more sophisticated variant. But I think the idea should be clear from the above.