According to physical intuition, it seems that the Sobolev norm shouldn't be a physically useful quantity. Why is this the case, and why isn't a more dimensionally correct mixed norm employed more often?
Consider a function $f: \mathbf{R}^d \to \mathbf{R}$. If we denote the spatial units of $f$ by $A$, and the units of the output of $f$ by $B$, then an easy calculation shows that the units of $\| f \|_p$ is equal to $A B^{d/p}$. On the other hand, the output of the $i$th partial derivative $f_i$ has units $A/B$, which means that $\|f_i\|_q$ has units $AB^{d/q-1}$. From the perspective of a dimensional analysis, $\| f \|_p + \| f_i \|_q$ should only make sense if $AB^{d/p} = AB^{d/q-1}$.
Why, in this case, do we not consider the mixed norm Sobolev spaces with $d/p = d/q-1$ more often than the standard Sobolev spaces $W^{1,p}(\mathbf{R}^d)$? Shouldn't the scale invariance make the mixed norm spaces more useful than the standard Sobolev spaces?
So first, quantities are not always quantities with units! For example if you want to speek about a Gaussian $e^{-|x|^2} = \sum_{n=0}^∞ \frac{(-1)^n\,|x|^{2n}}{n!}$, it clearly has no homogeneity. So, it means that the $x$ has to be adimensional. Usually, in physical models, you will get something like $e^{-λ\,|x|^2}$ where if the unit of $x$ is $L$, then $\lambda$ is of unit $L^{-2}$. In partial differential equations, a lot of equations are about a density of something, which often have no units. At the end, it will depend a lot on your specific physical system behind, and the questions you are asking about it.
About Sobolev norms, you could as well just multiply by a factor and get $λ\,\|f\|_{L^p} + \|∇f\|_{L ^p}$. In part of the papers, one prefers to take $\lambda=1$ (of course if you have units it could be $λ = 1$ m/s for example, so the dimensional analysis is lost (but easy to modify usually).
Then, as you say, these Sobolev spaces are not homogeneous. One usually defines homogeneous Sobolev spaces $\dot{W}^{1,q}$ as the completion of $C^{\infty}_c$ with respect to the norm $\|∇f\|_{L^q}$, and defines the norm $\|f\|_{\dot{W}^{1,q}} = \|∇f\|_{L^q}$ on this space. This is homogeneous, so you can put directly units. Most of the other spaces of functional analysis have a homogeneous version.
Remark actually that your dimensional analysis is a good way to find the exponents for the Sobolev embeddings, since $W^{1,q} ⊂ L^p$ with your exponents. In particular, it implies that the norm $\|f\|_{L^p} + \|∇f\|_{L ^q}$ is equivalent to the homogeneous Sobolev space $\dot{W}^{1,q}$, since $$ \|∇f\|_{L^q} ≤\|f\|_{L^p} + \|∇f\|_{L^q} ≤ (C+1)\,\|∇f\|_{L^q} $$ so there is no need to add the $L^q$ norm, unless one really is interested by the value of the constant, or some other precise effects.