# Numerical Solution of PDE with Singularity

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?

• How do you approximate $\frac{\partial P}{\partial x}$? I suppose the main reason for negative values are oscillations, if this is true, than you need to use a monotone approximation (e.g., upwind) for this term. – VorKir Jul 13 '17 at 16:45