Epsilon regularity: what does it say and where does it come from? The $\varepsilon$-regularity phenomenon shows up in several different contexts. I try to describe it focussing on the harmonic map situation, but I really would like to understand the situation in general.
The following is the Schoen-Uhlenbeck $\varepsilon$-regularity lemma, extracted from
Tobias H. Colding, William P. Minicozzi II, An excursion into geometric analysis.
Let $N$ be a Riemannian manifold and $B_{r}$ be the ball of radius $r$ centred at the origin in $\mathbf{R}^k$. Then there exists $\varepsilon(k,N)$ such that if $u:B_{r}\in\mathbf{R}^k\rightarrow N$ is an energy minimizing map and
$$\frac{\int_{B_{r}}|\nabla u|^2}{r^{k-2}}<\varepsilon,$$
then $u$ is smooth in a neighborhood of $0$ and
$$|\nabla u|^2(0)\leq \frac{C}{r}.$$
Thus if a (conformally invariant) rescaling of the energy that $u$ minimizes is small (I suppose $u$ should be in a suitable Sobolev space), then $u$ is automatically smooth in some smaller ball. This rescaling is monotonically increasing thanks to a monotonicity lemma. I am not sure how to interpret the bound on the derivative at zero, though.
The $\varepsilon$-regularity lemma quickly implies that the singular set $S$ of $u$ has $(k-2)$-dimensional Hausdorff measure zero.
My questions are:


*

* What are the basic ingredients (I suppose I am talking about the properties of the energy functional here) that guarantee that such a lemma holds?

* What is the meaning of the supremum of the set of all $\varepsilon$ such that the energy bound holds, and how can it be computed?

* Is there a simple intuitive picture that I am missing that explains the situation?

* Is there an instance of this phenomenon that predates the Schoen-Uhlenbeck paper?


Many thanks.
 A: I can comment on the $\epsilon$-regularity lemma for 4-dimensional Einstein manifolds. Namely, there is an $\epsilon$ depending on the dimension and Sobolev constant so that $\int_{B_x(r)} |Rm|^2 dV_g < \epsilon$ implies that $\sup_{B_x(r / 2)} |Rm|^2 \leq Cr^{-2} \int_{B_x(r)} |Rm|^2 dV_g$.
The key ingredient that makes this lemma work is that for Einstein manifolds, the function $|Rm|$ satisfies an elliptic inequality (it's ``subharmonic" for some elliptic operator). From there a standard PDE argument using Moser iteration gives the $\epsilon$-regularity. (It's like a non-linear version of the mean value inequality for subharmonic functions in Euclidean space).
A: Below is a rather longwinded description of the special case when the singularity is at worst an isolated point. I suspect you know all this already. The magic comes at the very end (see paragraph that starts with "Here's the critical trick"). I don't know if this is the same thing that gives Schoen-Uhlenbeck the extra oomph or not.
There are three applications I know of: Minimal hypersurfaces,
self-dual Yang-Mills connections, and Einstein manifolds. The
regularity theory described below is used for both a convergence
theorem except possibly a finite number of points and a removable
singularity theorem. These theorems are then used to establish the
so-called bubbling phenomenon. The story below applies to the latter
two applications; the details for minimal hypersurfaces are slightly different.
Assume for convenience that we're
on a smooth $n$-dimensional complete Riemannian manifold, where $n >
2$. Denote the
Laplacian on both functions and tensors by $\Delta = g^{ij}\nabla_i\nabla_j$.
Denote the $L_p$ norm of a function or tensor $u$ with respect to the
Riemannian metric by $\|u\|_p$.
Throughout the discussion below we will restrict to a geodesic ball
$B(x, r)$
and assume that the following Sobolev inequality holds for
a fixed constant $C_S$ and any smooth function $u$ compactly supported
in $B$:
$$
\|\nabla u\|_2 \ge C_S\|u\|_{2n/(n-2)}.
$$
First, you consider the scalar elliptic inequality $-\Delta u \le bu$,
where $b$ can be viewed as a given potential function. Using Moser
iteration, you show that if
$$ \|b\|_{q/2}, \|u\|_p < C, $$
where $q > n$, for
some $p > 1$ on $B(x,r)$, then there is a bound on
$\|u\|_\infty$ on, say, $B(x,r/2)$.
Second, you use Moser iteration to show that if $ \|b\|_{n/2} $ is sufficiently small (depends
on $C_S$) on $B(x,r)$, then there is a bound for $\|u\|_{q/2}$ for
some $q > n$ on $B(x,r/2)$.
Combining the first two shows that if $u$ satisfies $-\Delta u \le
cu^2$ and 
$ \|u\|_ {n/2} $
 is sufficiently small on $ B(x,r) $, then there
is a bound on $ \|u\|_\infty $ on $ B(x,r/2) $.
In each application there is a curvature tensor $F$ that satisfies a
PDE of the form
$$
-\Delta F = Q(F),
$$
where $Q$ depends quadratically on $F$. Moreover, there is a
convergence theorem when there is a uniform pointwise bound on $ F $
(for Einstein manifolds you use the Cheeger-Gromov convergence theorem).
Applying the results above to $u = |F|$ using coverings with smaller
and smaller balls leads to a convergence theorem when there is a
uniform bound on 
$ \|F\|_ {n/2} $ where the convergence can fail at only
a finite number of points (where in the limit there is too much of $ \|F\|_{n/2} $
for the estimates above to hold).
Now you want to study the limit object near each point
singularity. If you keep close track of the dependence on $r$ in the
estimates above, the best you can do is a bound on $F$ that blows up
like $r^{-2}$, where $r$ is the distance to the singularity. This is
not enough to remove the singularity, so you need to use more than the
elliptic PDE above.
Here's the critical trick: When doing Moser iteration on $u = |F|$,
you use the standard Cauchy-Schwarz inequality to obtain the following
pointwise inequality:
$$
|F\cdot\nabla F| \le |F||\nabla F|
$$
But in all of the applications, you have extra information about $F$ and its
covariant derivative. In particular,
$F$ and/or its covariant derivative have certain symmetries, which
allow you to prove a pointwise bound of the form
$$
|F\cdot\nabla F| \le c|F||\nabla F|,
$$
where $c < 1$. This improvement when used with Moser iteration allows
you to show that $F$ blows up more slowly than $r^{-2}$. Iterating
this improvement leads to a uniform pointwise bound on $F$, which in turn
allows the singularity to be removed using a straightforward geometric
ODE argument.
The removable singularity theorem allows you to analyze both the
limiting object with the bubbles removed as well as the bubbles themselves.
ADDED: I can't resist adding an anecdote to this: Right after I learned the trick in the paragraph above from a paper of Schoen-Simon-Yau, I went to a colleague's office to show it to him. As it happens, Eli Stein was there, and he exclaimed, "But it's in my book!" And indeed it is. You will find it presented very nicely in VII.3.1 "A subharmonic property of the gradient" of Stein's 1970 book, "Singular Integrals and Differentiability Properties of Functions". It is obvious that S-S-Y did not know this or forgot, because their proof is much messier than Stein's.
A: The way I think of it is to view semilinear PDEs, such as the harmonic map equation, as a contest between the linear portion of the equation ($\Delta u$ in this case) and the nonlinear portions (which, in the case of harmonic maps, are roughly of the shape $|\nabla u|^2$).  Intuitively, if the nonlinear part is small compared to the linear part then we expect the linear behaviour to dominate.  In the case of harmonic maps, this means that we expect the solutions to behave like solutions to Laplace's equation $\Delta u = 0$, which are known to be regular.
A bit of dimensional analysis then tells us that the condition $\frac{\int_{B_r} |\nabla u|^2}{r^{k-2}} < \varepsilon$ has the right scale-invariance properties to have a chance of making the nonlinear term smaller than the linear term.  (To make this rigorous, one of course needs to deploy various harmonic analysis estimates in well-chosen function space norms, such as Sobolev embedding.)
I discuss these heuristics (though more for dispersive equations than for elliptic ones) a bit at
http://terrytao.wordpress.com/2010/04/02/amplitude-frequency-dynamics-for-semilinear-dispersive-equations/
The question of what happens at the critical value of epsilon is an interesting one.  Often, the limiting non-regular solutions at that value of epsilon, after rescaling and taking limits, tend to be quite symmetric and smooth, away from a very simple singular set (e.g. a subspace).  I don't know the elliptic case too well, but one obvious candidate for such a solution would be a singular 2D harmonic map (such as the map from C -> S^1 given by x -> x/|x|) extended to k dimensions by adding k-2 dummy variables.  In the dispersive case, the analogous concept is that of the minimal energy blowup solutions, and these tend to be soliton solutions (so, typically, they obey a time translation invariance symmetry), associated to the ground state solution of the associated time-independent equation.
