Solution of Heat equation with Neumann BC in an arbitrary domain Consider the heat equation $u_t=\Delta u$ with Neumann boundary condition and initial condition $u(x,0)=u^0(x)$ in a bounded domain $\Omega$ with smooth boundary. 
Is this true:
Any solution $u(x,t)\in W^{2,p}$ of the equation can be written as $$u(x,t)=k(x,t)\star u^0(x)$$ where $k$ is a green function (depends on $\Omega$).  
 A: There exists a theory of Green functions for general parabolic boundary value problems which covers the case you are interested in, in particular papers by Eidelman, Ivasishen, Solonnikov. For references see 
S. D. Eidelman and N. V. Zhitarashu, Parabolic boundary value problems. Basel: Birkhäuser (1998).
Unfortunately, most of the papers on this subject are available only in Russian.
A: Probably yes, if I interpret your question correctly. The answer is based on the theory of operator semigroups.
First of all, let us consider the case $p=2$: Then the solution is given by a semigroup on $L^2(\Omega)$ generated by an operator associated with a quadratic form: it is known that the semigroup is ultracontractive (and in particular has a kernel $k$ that is of class $L^\infty(\Omega\times \Omega)$ with respect to space) if and only if the form domain $H^1(\Omega)$ is continuously embedded in some $L^p(\Omega)$ space for $p$ large enough, which luckily does hold by virtue of the Sobolev embedding theorem.
