Nice way of thinking about the Laplace operator... but what's the proof? Here's a nice fact: roughly speaking, the Laplace operator gives you the difference between the value of a function at a point and the average value at "neighboring" points.
More precisely, in $\mathbb{R}^n$ let $S_r(x)$ be the $(n-1)$-dimensional sphere of radius $r$ centered at the point $x$, let $V_r$ be the volume of this sphere, and let $d\sigma$ be the volume element on this sphere.  Then at every point $x \in \mathbb{R}^n$
$$ \lim_{r \rightarrow 0} \frac{\int_{S_r(x)} f d\sigma}{V_r} - f(x) = \frac{r^2}{2n} \Delta f(x) + \bar{o}(r^2) $$
for all $C^2$ functions $f$ on $\mathbb{R}^2$.  So far, however, I've been unable to prove this fact.  In the book I'm following (Grigor'yan, Heat Kernel and Analysis on Manifolds) the only theorem that's really been introduced so far is one of Green's identities:
$$ \langle u, \Delta v \rangle = \langle \nabla u, \nabla v \rangle, $$
where at least one of $u,v$ is compactly supported and $\langle \cdot, \cdot \rangle$ denotes the inner product over $\mathbb{R}^n$.  So, it seemed natural to consider any compactly-supported test function $g \in C^1(\mathbb{R}^n)$, in which case the formula above would look something like
$$ \frac{\int_{S_r(x)}\langle f, g \rangle d\sigma}{V_r} - \langle f, g \rangle = \frac{r^2}{2n} \langle \nabla f, \nabla g \rangle + \langle g + \bar{o}(r^2) \rangle. $$
Not sure where to take it from here, though (or even if this is the right direction!).  Any hints/tricks are much appreciated.  (For the record, I am not solving this problem as a homework exercise.)
Finally, a more minor question: what does the bar signify in $\bar{o}(r^2)$?
Thanks!
 A: Given a point $p$ in the domain of the function $f$, consider the function which is the average value $A(r)$ of $f$ over the sphere $S_r$ of radius $r$ centred at $p$.  If $S^{n-1}$ is the unit sphere in $\mathbb R^n$ then you can write the function as:
$$A(r) = \int_{S^{n-1}} f(p+rx)dx$$
I'm integrating with respect to the standard content on the sphere, but rescaled so that the sphere has unit content. 
Differentiating $A(r)$ with respect to $r$ gives
$$\frac{dA}{dr} = \int_{S^{n-1}} x\cdot\nabla f(p+rx)dx$$
in the above equation $x \in S^{n-1}$ is a unit vector, and the $\cdot$ is dot product of vectors. 
The above can be interpreted as a flux integral, so Gauss's theorem about flux integrals applies giving
$$\frac{dA}{dr} = r\int_{D^n} \nabla^2 f(p+rx)dx$$
now $x$ is a variable for the unit ball $D^n$. 
So this tells you the first few terms of the Taylor expansion for the function $A$, in particular $A(0)=f(p)$, $A'(0)=0$, and so on with the next terms being some multiple of $\nabla^2f(p)$.  A quick calculation says it should be $\frac{\nabla^2(f(p))}{n}$ but perhaps I've been too sketchy. 
