# 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?

PLEASE NOTE:

1. 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.

2. I thought of using Fourier transform once, but that would help only if the coefficients are constants, which is not true in my case.

3. I have gone through these links, please do not suggest me any of these answers: (a), (b), (c), (d)

• 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
• @VorKir Please provide more details – Mitradip Das Jul 18 '17 at 7:41