I'm looking for a numerical method based on finite differences which would persue me to solve Fokker-Planck equation: $$ \frac{\partial p}{\partial t} = -\sum_{i=1}^n\frac{\partial}{\partial x_i}(a_i(x)p)+\sum_{i,j=1}^n\frac{\partial^2}{\partial x_i \partial x_j }(b_{ij}(x)p)$$ with initial condition $$ p(x,0)=p_0(x) $$ and border condition $$ \lim\limits_{|x|\rightarrow \infty}p =0, $$ where $p=p(x,t)$ - unknown function, $a_i(x), \, b_{ij}(x)-$ smooth enough coefficients, and $n=2$ for my system.
The problem is that to apply any finite difference scheme, say, alternating direction implicit method, one should replace the unbounded domain with a closed bounded one, so the problem is solvable in this case. But that modified problem differs from the problem stated initially, and the numerical solution converges to the solution of another problem. Is it possible to prove that the numerical solution of the problem with bounded domain converges to the solution of the problem given above when radius of the domain goes to infinity? Or are there other numerical methods which allow to work with unbounded region?