# Neumann initial-boundary value problem for the Heat Equation using Duhamel's formula

How can I write down the solution to the Neumann initial-boundary value problem for the Heat Equation using Duhamel's formula?

$$u_t - \Delta u = f(t,x) \;\text{ in } \Omega \subset \mathbb{R}^n \;\text{bounded}\\ u = u_0(x) \;\,\text{at} \; \,t=0\\ \frac{\partial u}{\partial \vec{\nu}}=0 \;\, \text{at} \;\, \partial\Omega$$

By Duhamel's Principle, the solution to the PDE should be $$u(\boldsymbol{x},t)=\int_\Omega w(\boldsymbol{x}-\boldsymbol{y},t)u_0(\boldsymbol{y})\;d\boldsymbol{y}+\int_0^t\int_\Omega w(\boldsymbol{x}-\boldsymbol{y},t-s)f(\boldsymbol{y},s)\;d\boldsymbol{y}\;ds$$ where $w$ is known as the Neumann heat kernel and solves the following PDE \begin{align} w_t(\boldsymbol{x},\boldsymbol{y},t)-\Delta w(\boldsymbol{x},\boldsymbol{y},t)&=0\qquad &&(\boldsymbol{y},t)\in\Omega\times(0,\infty)\\ \lim_{t\to0}w(\boldsymbol{x},\boldsymbol{y},t)&=\delta(\boldsymbol{x}-\boldsymbol{y})\qquad &&\boldsymbol{y}\in\Omega\\ (\partial_\boldsymbol{\nu}w)(\boldsymbol{x},\boldsymbol{y},t)&=0\qquad &&\boldsymbol{y}\in\partial\Omega,\;t>0 \end{align} One can easily check this solution by plugging $u$ into your original PDE.