Removable Singularities for Elliptic Equations The following fact is quite standard and does not have a very long proof:
$(\*)$ If $u$ is harmonic on $B_1(0)\setminus \{0\}$ and uniformly bounded, then $u$ in fact extends to a harmonic function on the whole ball.
Some googling reveals that such statements are in fact true for large classes of elliptic operators with much more general singularity sets and growth conditions. There appears to be a vast literature on the subject.
I  would like to know if the exact analogue of $(\*)$ has a "simple" proof for linear elliptic operators. 
"Simple" is slightly ambiguous but is meant to mean tools present in standard elliptic theory textbooks, i.e. Gilbarg and Trudinger or Jost.
 A: Here are the details behind my original interpretation of Deane Yang's comment:
By Schauder estimates, $\vert\vert \nabla u\vert\vert_{L^{\infty}} \leq C\vert\vert u\vert\vert_{L^{\infty}}$, so we also know that the gradient of $u$ is uniformly bounded. Next, we claim that $u$ is a weak solution of $Lu = 0$ on $B_1(0)$. We have $Lu = \sum_{i,j=1}^n \partial_j(a^{ij}(x)\partial_i u) + \sum_{i=1}^nb^i(x)\partial_i u + c(x)u$. To show that $u$ is a weak solution on $B_1(0)$, we need to show that for every $\varphi \in C_c^{\infty}(B_1(0))$, 
$\int_{B_1(0)} \sum_{i,j=1}^n a^{ij}(x)\partial_{i}u\partial_j\varphi + \sum_{i=1}^nb^i(x)(\partial_i u) \varphi+ c(x)u\varphi = 0$
Since all of the relevant terms are uniformly bounded, the left hand side is equal to
$\lim_{\epsilon\to 0}\int_{B_1(0)\setminus B_{\epsilon}(0)} \sum_{i,j=1}^n a^{ij}(x)\partial_{i}u\partial_j\varphi + \sum_{i=1}^nb^i(x)(\partial_i u) \varphi+ c(x)u\varphi$
After an integration by parts, this gives
$\lim_{\epsilon\to 0}\int_{B_1(0)\setminus B_{\epsilon}(0)}(Lu)\varphi + \int_{\partial B_{\epsilon}(0)}\sum_{i,j=1}^n a^{ij}(x)(\partial_{i}u)\varphi \nu_j = $
$\lim_{\epsilon\to 0}\int_{\partial B_{\epsilon}(0)}\sum_{i,j=1}^n a^{ij}(x)(\partial_{i}u)\varphi \nu_j$
Since everything involved is uniformly bounded, this goes to $0$. (EDIT: The dimension should be at least $2$ for this to work)
Now $L^2$ elliptic regularity implies that $u$ is smooth on the whole ball and we are done.
A: The proof I had in mind was actually simpler and appears to work only in dimension greater than 2. Here is a sketch for a second order self-adjoint elliptic operator $Pu = \partial_i(a^{ij}\partial_ju)$. Suppose $u$ satisfies $Pu = 0$ on $B\backslash\{0\}$. We want to show that $u$ is a weak solution on $B$. It suffices to show that if $\phi$ is a smooth compactly supported function on $B$, then
$$
\int \phi Pu = 0.
$$
Let $\chi$ be a smooth nonnegative compactly supported function that is equal to $1$ on a neighborhood of $0$ and, given $\epsilon > 0$, let $\chi_\epsilon(x) = \chi(x/\epsilon)$. Then
$$
\int \phi Pu = \int [\phi(1-\chi_\epsilon) + \phi\chi_\epsilon]Pu = \int \phi\chi_\epsilon Pu
= \int P(\phi\chi_\epsilon)u.
$$
If we now expand out $P(\phi\chi_\epsilon)$, we see that it is $O(\epsilon^{-2})$. Therefore, if $u \in L_\infty$, then
$$
\int P(\phi\chi_\epsilon)u = O(\epsilon^{n-2})
$$
Therefore, taking the limit $\epsilon \rightarrow 0$ shows that $u$ is a weak solution on $B$. Elliptic regularity tells you that it is in fact a strong smooth solution.
A: Check out
Harvey, Reese; Polking, John Removable singularities of solutions of linear partial differential equations. Acta Math. 125 1970 39–56.
