Let $\phi(x)$ be a free massless scalar field on $d$-dimesnional space-time with Euclidean metric. I am interested in the operator formalism, i.e. $\phi(x)$ is an operator satisfying $\Delta \phi=0$ where $\Delta=\sum_{i=1}^d \partial_{i}^2$.
I would like to compute the two point function $$D(x):=\langle 0|\phi(x)\phi(0)|0\rangle$$ (there is no time ordering!).
The following elementary argument below leads to a conclusion which does not sound to me to make sense, namely it proves that $D(x)$ is a constant function. My question is whether this answer is indeed correct.
Since $\Delta\phi(x)=0$ one has $\Delta D(x)=0$. Moreover $D(x)$ is rotation invariant, i.e. depends on $r$ only. The space of such $SO(d)$-invariant solutions outside of 0 must be a linear combination of the constant function 1 and either $\frac{1}{r^{d-2}}$ when $d>2$ or $\log (r)$ when $d=2$. Since $\Delta \frac{1}{r^{d-2}}$ is proportional to $\delta(x)$ (similarly in 2d) the only possibility is that $D(x)$ is the constant function. This last conclusion does not sound to make sense.
Remark. If one considers the same question in Minkowski space, say in 4d, the 2-point function is $$D(x)=\int \frac{d^4p}{(2\pi)^3}\delta(p^2)\theta(p_0)e^{-ip\cdot x}=\int\frac{d^3p}{(2\pi)^3}\frac{1}{2|\vec p|}e^{-i(|\vec p|t-\vec p\cdot \vec x)}.$$ This function does satisfy $\Box D(x)=0$ and is Lorentz invariant.