Approximation argument for Poisson integral formula In the PDE book,

for a harmonic function $u \in C^2(B_R(0))\cap C^1(\overline B_R(0)),$ we have the following Poisson integral formula
  $$
u(y)=\frac{R^2-|y|^2}{n w_n R}\int_{\partial B_R(y)}\frac{u(x)}{|x-y|^n}ds_x
$$
  An approximation argument shows that the Poisson integral formula continues to hold for $u \in C^2(B_R(0))\cap C(\overline B_R(0)).$

What is the approximation argument? For a harmonic function $u \in C(\overline B_R(0))$, is there a sequence of harmonic functions $u_k \in  C^1(\overline B_R(0))$ such that $u_k \to u$ in $C(\overline B_R(0))$? I have no idea for harmonic functions. May I use convolution?
Please let me know if you have any hint or comment for it. Thanks in advance!
 A: A typical way to approximate a function defined on a ball $B_R(0)$ is 
$$u_k(x) = u(\lambda_k x),\quad \lambda_k\nearrow 1$$
For any $u \in C(\overline B_R(0))$, this results in $u_k \to u$ in $C(\overline B_R(0))$ by virtue of the uniform continuity of $u$.  If $u$ is smooth (harmonic, holomorphic...) in the open ball $B_R(0)$, then $u_k$ has these properties  on a larger ball.
A: There is a sequence of polynomials converging to $u$ in $C(\overline {B}_R(y))$ by Stone Wierstrass Theorem.
A: For $t\in(0,1)$, consider $\lambda: \overline{B_r(0)}\rightarrow B_r(0)$, $\lambda(y)=ty$, then $\lambda\in C^{\infty}(\overline{B_r(0)})$.  As $u\in C^2(B_r(0))$ follows that $\hat{u}(y)=u(ty)\in C^2(\overline{B_r(0)})$. Moreover
$$
\frac{u(ty)}{|x-y|^n}\rightarrow \frac{u(y)}{|x-y|^n}, \text{ as } t\rightarrow 1, \text{ and } 
\frac{|u(ty)|}{|x-y|^n}
\leq \frac{\displaystyle\max_{\overline{B_r(0)}}|u|}{|x-y|^n}\in L^1(\partial B_r(0)).
$$
Therefore, from the dominated convergence theorem, for all $y\in B_r(0)$,
\begin{eqnarray*}
u(y)&=&\lim\limits_{t\rightarrow 1}u(ty)
=\lim\limits_{t\rightarrow 1}\hat{u}(y)\\
&=&\frac{r^2-|y|^2}{nw_nr} \lim\limits_{t\rightarrow 1}\int_{\partial B_r(0)}\frac{u(tx)}{|x-y|^n}ds_x
=\frac{r^2-|y|^2}{nw_nr} \int_{\partial B_r(0)}\lim\limits_{t\rightarrow 1}\frac{u(tx)}{|x-y|^n}ds_x\\
&=&\frac{r^2-|y|^2}{nw_nr} \int_{\partial B_r(0)}\frac{u(x)}{|x-y|^n}ds_x.
\end{eqnarray*}
