Given the heat problem:
$$\begin{cases} \frac{d}{dt}u(x,t)=\Delta u(x,t) & \forall (x,t)\in \Omega\times(0,T) \\ u(x,0)=u_0(x) & \forall x\in\Omega \\ u(x,t)=0 & \forall x\in\partial\Omega, t>0 \end{cases}.$$ if $\Omega$ is smooth enough (not necessarily bounded) and $u_0\in L^2(\Omega)$, one can assure that we have solutions $u\in C^\infty(\bar{\Omega}\times(0,T))$ (the initial data is reached in the sense that $u(t)\to u_0$ in $L^2(\Omega)$ when $t\to 0$.) This can be proved, for example, with Semigroup Theory: $(-\Delta,H^2(\Omega)\cap H^1_0(\Omega))$ generates a semigroup in $L^2(\Omega)$ ([Davies] Sections 1.3, 1.4) which is analytic. Therefore, $u(t)\in D((-\Delta)^{k})$ for every $k\in \mathbb{N}$. Now, using higher regularity results (see for example Evans Section 6.3 Theorem 5) one can prove that $D((-\Delta)^k)\subset H^{2k}(\Omega)$ for every $k\in \mathbb{N}$. Therefore, $u(t)\in C^\infty(\bar{\Omega})$ for every $t>0$. However, one can consider a more general problem: $$\begin{cases} \frac{d}{dt}u(x,t)=\Delta u(x,t) & \forall (x,t)\in \Omega\times(0,T) \\ u(x,0)=u_0(x) & \forall x\in\Omega \\ Bu(x,t)=0 \forall & x\in\partial\Omega, t>0 \end{cases}.$$ where $Bu$ corresponds to Dirichlet ($Bu=u$), Neumann ($Bu=\frac{\partial u}{\partial n}$) or Robin ($Bu=\alpha\frac{\partial u}{\partial n}+\beta u$) and $u_0\in L^p(\Omega)$ where $1\leq p \leq \infty$. I think that, if $\partial\Omega$ and the cofficients of $B$ are smooth, it is well-known that $u\in C^\infty(\bar{\Omega}\times(0,T))$. However, I do not find any good reference for this.
In this general case, I know that $-\Delta$ (with an adequate domain, see [Denk]) generates a semigroup in $L^p(\Omega)$ (See for example [Davies]). However, it is not analytic if $p=1,\infty$. Furthermore, there are regularity results in the $L^p$ framework (See [GT] Theorem 9.19. In this case higher regularity in $L^p$ is left as an exercise in the homogeneous Dirichlet case (I think the authors forgot to mention that it is for Dirichlet)). There are also regularity results for the Poisson problem with other boundary conditions (See for example [Mikhailov] P.217 Theorem 4, which is stated for Dirichlet and Neumann, and the footnote includes Robin).
However, I could not deduce from all these result the fact that $u\in C^\infty(\bar{\Omega}\times(0,T)$ in the general case (that is, when $u_0\in L^p(\Omega)$ and for Dirichlet, Robin or Neumann B.C.) but just for some particular cases as $1<p<\infty$ with Dirichlet boundary conditions or $p=2$ with any boundary condition. I feel there may be an easier way to do this. Do you know how? Or, do you know any reference for this?
References:
[Evans] Partial differential equations. Evans, L. C. (2022).
[Davies] Heat Kernels and Spectral Theory. E. B. Davies.
[Denk] New thoughts on old results of RT Seeley. R Denk, G Dore, M Hieber, J Prüss, A Venni.
[GT] Elliptic Partial Differential Equations of Second Order. D Gilbarg, NS Trudinger.
[Mikhailov] Partial Differential Equations. V.P. Mikhailov.