Suppose that $\rho \in L^1(\mathbb{R}^n \times (0,T))$ for every $T < \infty$ is a weak solution of the PDE \begin{align} \partial_t\rho &= \Delta \rho + \text{div}(\rho\nabla\Psi(x))\\ \rho(t = 0) &= \rho^0 \end{align} where $\Psi \in C^\infty(\mathbb{R}^n)$ is given and $\Delta$ denotes the Laplacian; that is, \begin{equation} \int_{0}^{\infty} \int_{\mathbb{R}^n} \rho\left(\partial\xi - \nabla\Psi \cdot\nabla\xi + \Delta\xi\right)d x d t + \int_{\mathbb{R}^n}\xi(0)\rho^0d x = 0 \end{equation} for all $\xi \in C_c^\infty(\mathbb{R}^n\times\mathbb{R})$.

I am trying to understand how to prove that this weak solution $\rho$ is in fact a classical smooth solution.

The paper "The Variational Formulation of the Fokker-Planck Equation" http://www.imati.cnr.it/savare/Ravello2010/JKO.pdf sketches this argument, but I am having trouble understanding one part. In particular, they obtain the following expression on page 15, Equation (55): \begin{align} (\rho\eta)(t_1) &= \int_{t_0}^{t_1}\left[\rho(t)(\Delta\eta - \nabla\Psi\cdot\nabla\eta)\right]*G(t_1-t)dt\\ &\quad + \int_{t_0}^{t_1}\left[\rho(t)(2\nabla\eta - \eta\nabla\Psi\right]*\nabla G(t_1 - 1)dt\\ &\quad + (\rho\eta)*G(t_1 - t_0) \end{align} for all $\eta \in C_c^\infty(\mathbb{R}^n)$ and for a.e. $0 \leq t_0 < t_1$, where $G(x,t) = \frac{1}{(2\pi t)^{n/2}}e^{-|x|^2/2t}$ is the heat kernel.

After a few straightforward computations, they show that $\rho \in L^p_{loc}(\mathbb{R}^n\times(0,\infty))$, $p < \frac{n}{n-1}$.

Then, in the following line, all they say is "We now appeal to the $L^p$-estimates [18, section 3, (3.1), and (3.2)] for the potentials in (55) - [the above integral equality] - to conclude by the usual bootstrap arguments that any derivative of $\rho$ is in $L^p_{loc}(\mathbb{R}^n\times(0,\infty))$, from which we obtain the stated regularity condition ($\rho \in C^\infty(\mathbb{R}^n\times(0,\infty)$).

Could someone provide a good explanation generally for how a bootstrapping argument would work with $L^p$ estimates, and specifically in this particular case. Also how do $L^p$ estimates from reference [18], (which is O. A. Ladyzenskaja, V. A. Solonnikov, and N. N. Uralceva, Linear and Quasi–Linear Equations of Parabolic Type, American Mathematical Society, Providence, RI, 1968) come into play.

Any elucidation on this or suggestions for reading would be greatly appreciated!

NOTE: This question was originally posted on Math StackExcahnge a few months ago (https://math.stackexchange.com/questions/843901/regularity-of-a-weak-solution) and a few days ago by me (https://math.stackexchange.com/questions/941505/regularity-of-the-fokker-planck-equation?lq=1) which I have deleted since I came to know of the earlier post. Since there have been no replies I have slightly modified the earlier question and posted it here.