How should I set the boundary condition of a diffuse equation in a potential well?

Question

I want to solve the initial value problem of a particle diffusing in a well on a disk.

Paticle density$\rho(\vec{r},t)$satisfy the following equation: $$ \partial\rho/\partial t=\nabla\cdot(\nabla V(r)\rho)+D_t\nabla_{x,y} ^2\rho\\ \frac{\partial\rho}{\partial\hat{n}}=0,\text{on the boundary } (|\vec{r}|==1)\\ \rho(r,0)=\text{Gaussian}(r) $$ $V(r)$ is a potential well, $V(r)=r^2$.

I try to solve this problem in Mathematica with the Finite element method.

But the integral of the particle density does not conserve (there appears flux on the boundary).

It seems to be due to my false boundary condition.

How should I set the boundary condition?

Answer 1

The diffusion equation you are using can written as a balance law, a simple continuity equation $$ \partial_t \rho + \nabla \cdot \vec{j} = 0\,, $$ using a current density $$ \vec{j} = - \nabla V \rho - D_t \nabla \rho\,. $$ For the mass inside a volume you integrate the differential balance law to obtain an integral balance law $$ \partial_t \int_\Omega \rho \text{d}V = - \int_\Omega \nabla \cdot \vec{j} \text{d}V = - \int_{\partial\Omega} \vec{j} \cdot \vec{n} \text{d}S\,. $$ The very basics of continuum mechanics... ;) In words: The rate of change of the mass/probability in $\Omega$ (left) equals the outflux (right) that is computed by integrating the current density normal to the surface.

To ensure conservation in your problem, you thus need to make sure that the integral on the right vanishes, for example by making $\vec{j}$ vanish everywhere on the boundary.