Shouldn't also the Sobolev norm exist (be finite)?
Yes, but this condition is automatically satisfied because of the way that the space was defined. In fact, by definition
$$L_2(\Omega)=\left \{f:\Omega\to\mathbb{R};\;f\text{ is measurable and }\int_\Omega |f(x)|^2 \ dx<\infty\right\}.$$
So, "$u \in L_2(\Omega)$" implies
$$\int_\Omega |u(x)|^2 \ dx<\infty$$
and "$\partial^\alpha u \in L_2(\Omega)$ for all $|\alpha|\leq k$" implies
$$\int_\Omega |\partial^\alpha u(x)|^2 \ dx<\infty,\quad \forall\ |\alpha|\leq k.$$
Thus, your definition implies that $\|u\|_{H^k(\Omega)}$ is finite (that is, the "existence" (finiteness) of the Sobolev norm is implicitly stated).
Remark: When you said and so its weak derivative, you were using the following property:
If $u$ possess classical derivative $u'$, then $u$ possess weak derivative and thus $u\in H^1(0,1)$. Moreover, $u'$ coincides with the weak derivative of $u$.
However, by definition of $H^1(0,1)$, the weak derivative of a function in $H^1(0,1)$ have to be a function in $L_2(0,1)$ and thus the above property is false because sometimes $u'$ exists but doesn't belong to $L_2(0,1)$ (that is, doesn't satisfy
$\int_0^1 |u'(x)|^2 \ dx<\infty.$
An example is given by your function $u(x)=x^{-1/4}$).
The valid property is:
If $u \in C^1(0,1)\cap L_2(0,1)$ and if $u'\in L_2(0,1)$, then $u\in H^1(0,1)$ . Moreover, $u'$ coincides with the weak derivative of $u$.
Just for convenience, let $K:=[a,b]^n$ for $a<b$ s.t. supp $\phi\subseteq K$. By the steps you already made in the question, we are left to bound
$$\| u\|_{p_i}^{p_i}=\int_K dy~ |u(y)|^{p_i}\,,$$
for any $i$. Let us first assume that $u\in C_c^\infty(\Omega)\subseteq C_c^\infty(\mathbb{R}^n)$. With the neat shorthand notation $((x,y)):=(y_1,\ldots,y_{i-1},x,y_i,\ldots,y_{n-1}))$, we have
\begin{align}
\int_K dy~ |u(y)|^{p_i}&=\int_a^b dx\int_{[a,b]^{n-1}}dy|u((x,y))|^{p_i}\\
&\leq(b-a)\sup_{x'\in[a,b]}\int_{[a,b]^{n-1}}dy~|u((x',y))|^{p_i}\\
&\leq(b-a)\sup_{x'\in[a,b]}\int_a^{x'} dx\int_{[a,b]^{n-1}}dy~ p_i|u((x,y))|^{p_i-1}|\partial_i u((x,y))|\\
&\leq (b-a)p_i\int_K dy~|u(y)|^{p_i-1}|\partial_i u(y)|\,,
\end{align}
where, in the second to last step, we have used the fundamental theorem of calculus, together with the fact that $\partial_i|f|\leq|\partial_if|$ where $f$ stays away from zero.
Let us denote $\| f\|_{p,K}=(\int_K |f|^p)^{1/p}$. When we use Hölder's inequality on the above, we find
$$\| u\|_{p_i,K}^{p_i}\leq \left\||u|^{p_i-1}\right\|_{q,K}\Vert\partial_iu\|_{p_i,K}$$
for $q$ with $\frac{1}{q}+\frac{1}{p_i}=1$. Since this implies $q(p_i-1)=p_i$ as well as $p_i/q=p_i-1$, the definitions give us
$$\| |u|^{p_i-1}\|_{q,K}=\|u\|_{p_i,K}^{p_i-1}\,.$$
The following trick is a true classic! Note that we have obtained
$$\|u\|^{p_i}_{p_i,K}\leq (b-a)p_i\|u\|_{p_i,K}^{p_i-1}\|\partial_iu\|_{p_i},$$
in which we can move $\| u\|_{p_i,K}^{p_i-1}$ to the left, and find
$$\|u\|_{p_i,K}\leq (b-a)p_i\|\partial_iu\|_{p_i}\,.$$
Remember, we only proved this inequality for $u\in C_c^\infty(\Omega)$. Still, it implies that a Cauchy sequence $(u_n)$ in $D^{1,\vec{p}}$ is a Cauchy sequence w.r.t. $\|\cdot\|_{p_i,K}$ as well, implying that $\|\lim u_n\|_{p_i,K}<\infty$. Thus we have $\|u\|_{p_i,K}<\infty$ for any $u\in D^{1,\vec{p}}$, which means the answer to your question is... yes.
Best Answer
Guessing
First I would make a quick guess based on the chart of function space. It groups together Sobolev spaces $W^{s,p}$ with the same value of $\frac n p -s$, because these are related by the embedding theorem. While the inclusion provided by this theorem is strict, the sharpness of the theorem still makes "the spaces with equal $\frac n p -s$ are similar" a useful heuristic.
Your function $f$ fails to be in $W^{1,1}$, because its gradient is not an $L^1$ function but rather a vector-valued measure supported on the boundary of $A$. On the other hand, this measure has finite mass (meaning $f\in BV$), which is pretty close to $W^{1,1}$. So, it seems that $s=1$, $p=1$ is at the edge of spaces to which $f$ belongs. From $$\frac{n}{1}-1 = \frac{n}{2} - s$$ we conclude that for $p=2$ (your question), the relevant $s$ is $s=1-\frac n2$.
This agrees with the situation in 1D, and suggests we won't find anything good when $n\ge 2$. Of course, this is just a guess; it may well be wrong.
Estimating
Recall the layercake principle: the integral of a nonnegative function $g$ is equal to $\int_0^\infty |\{g>t\}|\,dt$ where $|\cdot|$ stands for the measure of the set. So, let's consider the inequality $$ \frac{|f(x)-f(y)|}{\|x-y\|^{2s+n}}>t \tag{1} $$ for large values of $t$ (only they are of interest when the measure space has finite measure).
For (1) to hold, exactly one of $x,y$ must be in $A$; also, both must be within distance $\delta \approx t^{-1/( 2s+n)}$ of the boundary of $A$. This constrains $x$ to a set of measure $\approx \delta $. Also, $y$ must lie in a ball of radius $\delta$ around $x$. So,
$$ \left|\left\{ (x,y) : \frac{|f(x)-f(y)}{\|x-y\|^{2s+n}}>t\right\}\right| \approx \delta^2 \approx t^{-(n+1)/(2s+n)} $$ The integral over $t\ge 1$ converges iff $$\frac{ n+1}{2s+n} >1$$ which is equivalent to $s<1/2$.
(In particular, the "guess" wasn't correct after all.)