I am trying to solve a PDE in $x$ and $t$, which my equation looks as
$$\frac{\partial P}{\partial t}=A(x,t)\frac{\partial^2 P}{\partial x^2}+B(x,t)\frac{\partial P}{\partial x}+C(x,t)P+D(x,t)$$
where $P(x,t)$ represents a probability distribution (I am trying to solve a modified version of Smoluchowski equation). I used Crank Nicolson Algorithm, but since my starting point is a Dirac Delta function (estimated with very thin Gaussian/Lorentzian), I am getting a negative value of $P(x,t)$ at some points.
With the doubt that there may be error in my code, I calculated the values by using excel by employing the Euler method and found that the problem lies in the almost-singular Dirac Delta function (estimated with very thin Gaussian/Lorentzian). Is there some method to handle this?
PLEASE NOTE:
Please do not suggest me to go for analytical solution, since the equation I am working on is not solvable (or rather, has not been solved till now) analytically.
I thought of using Fourier transform once, but that would help only if the coefficients are constants, which is not true in my case.
I have gone through these links, please do not suggest me any of these answers: (a), (b), (c), (d)
