I don't understand the different behaviour of the advection-diffusion equation when I apply different boundary conditions. My motivation is the simulation of a real physical quantity (particle density) under diffusion and advection. Particle density should be conserved in the interior unless it flows out from the edges. By this logic, if I enforce Neumann boundary conditions the ends of the system such as $\frac{\partial u}{\partial x}=0$ (on the left and the right sides) then the system should be "closed" i.e. if the flux at the boundary is zero then no particles can escape.
For all the simulations below, I have applied the Crank-Nicolson discretization to the advection-diffusion equation and all simulation have $\frac{\partial u}{\partial x}=0$ boundary conditions. However, for the first and last rows of the matrix (the boundary condition rows) I allow $\beta$ to be changed independently of the interior value. This allows the end points to be fully implicit.
Below I discuss 4 different configurations, only one of them is what I expected. At the end I discuss my implementation.
Diffusion only limit
Here the advection terms are turned off by setting the velocity to zero.
Diffusion only, with $\boldsymbol{\beta}$=0.5 (Crank-Niscolson) at all points
The quantity is not conserved as can be seen by the pulse area reducing.
Diffusion only, with $\boldsymbol{\beta}$=0.5 (Crank-Niscolson) at interior points, and $\boldsymbol{\beta}$=1 (full implicit) at the boundaries
By using fully implicit equation on the boundaries I achieve what I expect: no particles escape. You can see this by the area being conserved as the particle diffuse. Why should the choice of $\beta$ at the boundary points influence the physics of the situation? Is this a bug or expected?
Diffusion and advection
When the advection term is included, the value of $\beta$ at the boundaries does not seem to influence the solution. However, for all cases when the boundaries seem to be "open" i.e. particles can escape the boundaries. Why is this the case? Is it because, as pointed out here, that the advection equation in the form I have written it doesn't allow propagation in the reverse direction so the wave just passes through?
Advection and Diffusion with $\boldsymbol{\beta}$=0.5 (Crank-Niscolson) at all points
Advection and Diffusion with $\boldsymbol{\beta}$=0.5 (Crank-Niscolson) at interior points, and $\boldsymbol{\beta}$=1 (full implicit) at the boundaries
Implementation to follow, but it is this diffusion term implementation, combined with this advection term.