Problem
I have to solve the nonhomogeneous classic problem
$$\left(P_{1}\right)\;\,\left\{ \begin{aligned} u_{t}\;-\; \Delta u\; &= \;f& &\textrm{on}\;\;\; \Omega\times\left(0,\,\infty\right) &\\ u\; &=\; 0 & &\textrm{on}\;\;\; \partial\Omega\times\left(0,\,\infty\right)&\\ u\; &= \;g & &\textrm{on}\;\;\; \Omega\times\left\{t\;=\;0\right\}& \end{aligned}\right.\\$$
where $\,\Omega = B(x_0,R).\,$ We assume that $\,f=f\left(x,t\right)\;$ and $\,g=g\left(x,t\right)\,$ are smooth enough.
My attempt
Following Lawrence C. Evans guidelines (page 45), I first tried to solve the homogeneous equation looking for a particular solution:
$$u\left(x,t\right)=\frac{1}{t^\alpha}v\left(y\right),\;y=\frac{x}{t^\beta}$$
Once I assumed the solution to be radial, $\,v\left(y\right)=w\left(\left|y\right|\right)=w\left(r\right)\;$ I solved the ODE:
$$\frac{w}{2}r^n\,+\,w'r^{n-1} = A$$ $$\displaystyle w\left(r\right)=Be^{-\frac{r^2}{4}}\,+\,Ae^{-\frac{r^2}{4}}\cdot\int_{0}^{r}\frac{e^{\frac{s^2}{4}}}{s^{n-1}}ds$$
I got stuck just here. I cannot assume $\,A=0\,$ and obtain the fundamental solution $\,\Phi\,$ because it doesn't verify the BC's.
How can I go on? Any kind of help would be useful. Thanks in advance!
The spherical symmetric solution in your question provides a basis for the PDE in the infinite space $R^n$. Your boundary value problem in a sphere of a finite radius, you need to expand the solution in the series of products of the spherical Bessel function of $r$ and $n$-dimensional spherical harmonics.
The Laplace operator in $3$-dimensional spherical coordinate is
$$0=-\frac{\partial u}{\partial t}+\Delta u =\frac{\partial u}{\partial t}-\frac{1}{r^2}\frac{\partial}{\partial r} \left(r^2 \frac{\partial u}{\partial r}\right) + \frac{1}{r^2 \sin\theta} \frac{\partial}{\partial \theta} \left(\sin\theta \frac{\partial u}{\partial \theta}\right) + \frac{1}{r^2 \sin^2\theta} \frac{\partial^2 u}{\partial \varphi^2} \tag1$$ Suppose $u(r,\theta,\phi)=T(t)R(r)\Theta(\theta)\Phi(\varphi)$. Substitute it into Equation (1) and divide it by $f$, we have
$$0=-\frac1T\frac{\partial T}{\partial t}+\frac{1}{r^2R}\frac{\partial}{\partial r} \left(r^2 \frac{\partial R}{\partial r}\right) + \frac{1}{r^2\Theta \sin\theta} \frac{\partial}{\partial \theta} \left(\sin\theta \frac{\partial \Theta}{\partial \theta}\right) + \frac{1}{r^2 \sin^2\theta\Phi} \frac{\partial^2 \Phi}{\partial \varphi^2} \tag2$$ We see that the first term is independent of $(r,\theta,\varphi)$ while the remaining terms are independent of $t$. The first term must then be a constant. Let us call it $-\lambda$. We have $$T'(t)=-\lambda T(t) \tag3$$ Multiplying (2) by $r^2$, we see that the third and fourth terms are independent of $r$ while the second term depends only on $r$. So we introduce another constant $-\mu$. $$\frac{1}{r^2R}\frac{d}{d r} \left(r^2 \frac{d R}{d r}\right)+\Big(\lambda-\frac\mu {r^2}\Big)R=0. \tag4$$ Apply the same trick again, we obtain $$\Phi''(\varphi)+\nu\Phi(\varphi)=0 \tag5$$ for some constant $\nu$ and for $\Theta(\theta)$ $$\Theta''+\cot\theta\,\Phi'+(\mu-\nu(\csc\theta)^2)\Theta=0 \tag6$$ for some constant $\mu$.
We will solve Equation ($4$) for with boundary condition $R(r=1)=0$ where we set the radius of the sphere at $1$.
...
Jumping to the conclusion first, let $$q_{\lambda,k,m}(r,\theta,\varphi):=j_k(r\sqrt\lambda)Y_{k,m}(\theta,\varphi),\ \lambda>0,$$ where $j_k$ is the spherical Bessel function with $\sqrt\lambda$ being a positive root of it, $Y_{k,m}$ is the spherical harmonic function. Expand $$g(r,\theta,\varphi)=\sum_{\lambda,|m|\le k}a_{\lambda,k,m}q_{\lambda,k,m}(r,\theta,\varphi)$$ to solve for the coefficient $a_{\lambda,k,m}$. This is facilitated by the property that the set $\{q_{\lambda,k,m}(r,\theta,\varphi)\}_{\lambda,k,m}$ is orthogonal in $L^2$.
In the case of the homogeneous PDE where $f=0$, $$u(r,\theta,\varphi,t) = \sum_{\lambda,|m|\le k}a_{\lambda,k,m}q_{\lambda,k,m}(r,\theta,\varphi)e^{-\lambda t}.$$ By Duhamel's principle, the inhomogeneous equation can be viewed as linear superposition of the homogeneous problem and $$u(r,\theta,\varphi,t) = \sum_{\lambda,|m|\le k}a_{\lambda,k,m}q_{\lambda,k,m}(r,\theta,\varphi,t)e^{-\lambda t}+\int_0^t ds \sum_{\lambda,|m|\le k}a_{\lambda,k,m}(s)q_{\lambda,k,m}(r,\theta,\varphi)e^{-\lambda (t-s)},$$ where $$f(r,\theta,\varphi,t)=\sum_{\lambda,|m|\le k}a_{\lambda,k,m}(t)q_{\lambda,k,m}(r,\theta,\varphi).$$
This is how one builds up the solution in a sphere. More details can be found in e.g. ON GREEN'S FUNCTIONS IN THE THEORY OF HEAT CONDUCTION IN SPHERICAL COORDINATES by Arnold N. Lowan and Green’s Function for the Heat Equation by Abdelgabar Adam Hassan.
Now if $g(x)$ is spherically symmetric thus is only a radial function $g(r)$, then the solution is simpler.