 Solution

The interior elliptic regularity assertion for the Laplace equation is that a harmonic function is smooth, in fact real analytic, throughout $U$. No boundary regularity of its unspecified boundary values is implied.

First let $u\in C^2(U)$ and $\Delta u=0$. On a ball compactly contained in $U$, differentiating its spherical average and applying the divergence theorem expresses that derivative as a constant factor times ${\int }_{B_r}\Delta u$, which is zero. The spherical average tends to $u$ at the center as $r\to 0$. This proves the mean value property for harmonic functions.

Choose a radially symmetric smooth mollifier ${\eta }_r$, supported in $B_r$, with integral one. By integrating the spherical mean value property,whenever $B_r(x)⋐U$. For fixed $r$ the right side is a smooth convolution, since all derivatives can be placed on ${\eta }_r$. Thus $u$ is smooth. The same argument proves the Weyl lemma for a distributionally harmonic $u\in L_{\mathrm{loc}}^1$: first mollify $u$, apply the fixed-radius identity, and let the mollification radius tend to zero in distributions to obtain the same smooth representative.

To prove real analytic regularity, differentiating the fixed-radius convolution gives the interior derivative estimate for a harmonic functionfor any harmonic $v$. All derivatives of $u$ are harmonic. On nested balls between $B_R(a)$ and $B_{R/2}(a)$, apply this estimate $k$ times, decreasing the radius by $R/(2k)$ each time. For $|\alpha |=k$,The inequality $k^k\le e^{k}k!$ follows by integrating $\log s$ below the sum defining $\log (k!)$. Apply the one-dimensional Taylor theorem along each segment, expanding directional derivatives by the multinomial formula. The remainder is bounded by $M(A{\sum }_j|h_j|{)}^k$, so it tends to zero for sufficiently small $h$. This gives a locally convergent multivariate Taylor series. A harmonic function is real analytic in the interior.

The usual inhomogeneous elliptic regularity statement also follows: if $\Delta u=f$ and $f$ is smooth, take a cutoff $\chi$ equal to one near a given point and set $w=\Phi \ast (\chi f)$, where $\Phi$ is a fundamental solution of the Laplace equation with $\Delta \Phi =\delta$. Moving every derivative to the compactly supported smooth function $\chi f$ shows $w$ is smooth. Locally $u-w$ is harmonic, so $u$ is smooth there too.