I read this post, "Conservation of a physical quantity when using Neumann boundary conditions applied to the advection-diffusion equation" and although it is the same type of equation it does not fit my problem.

I try to solve the following equation with a no flux boundary with the solver d03pf from NAG. Should not be difficult, but i do have troubles with a no flux boundary.

$ \frac{\partial(u)}{\partial(t)}+\frac{\partial(F(u))}{\partial(x)}=\frac{\partial}{\partial(x)}\left( D(u)\frac{\partial(u)}{\partial(x)}\right)$

I do have a closed vessel and particles are sedimenting towards the bottom. Before they touch each other it is of hyperbolic type. Otherwise the PDE is parabolic.

My flux function is something like $F=u\times\text{velocity}(1-u)^C$ So when i use the standard no flux function for the hyperbolic part, saying the flux at the boundary is zero. The solver sets the mass in the last cell to zero and i am loosing mass.

I can say $\text{flux}(\text{end}-1) = -\text{flux}(\text{end})$ and the solver raises the concentration in the last cell greater one (for e.g C=5). But that does not work for even C.

$F(u)$ as well as $D(u)$ are piecewise differentiable.

I do have the space vector, the time at this iteration and the solution of the last three cells given.

Does anybody know another way constructing a no flux boundary ?

Best wishes


  • $\begingroup$ Does your solver supports Robin condition? $\endgroup$
    – Dr_Sam
    Commented Mar 21, 2013 at 11:47
  • $\begingroup$ Hey Moritz is it possible that you can maybe provide us with the code of what you did to solve the issue? Thanks $\endgroup$
    – user17923
    Commented Oct 18, 2015 at 6:17
  • $\begingroup$ It is a long time ago but yes I can (if I find my backup folders). Give me a couple of days. $\endgroup$
    – Moritz
    Commented Oct 18, 2015 at 7:28
  • $\begingroup$ Welcome to SciComp.SE. I hope Moritz is able to help you but this should really be a comment and not an answer. $\endgroup$
    – James
    Commented Oct 18, 2015 at 7:40
  • $\begingroup$ Thanks a ton Moritz...I really appreciate it, Ive been reading all your threads on this forum and on mathworks. I've been having the same issues as you, looking forward to your reply! $\endgroup$
    – user17923
    Commented Oct 18, 2015 at 17:18

3 Answers 3


The correct boundary condition that you need to setup is

$$ F(u) - D(u) \frac{\partial u}{\partial x} = 0$$

If you could have this condition, the mass would be perfectly conserved (see the answer to the question that you refer to for the reasoning). The problem is that what this condition is non-linear, so this is not a Robin condition.

Here is what I would try to overcome this issue: make one part of the condition explicit in time. One option is to use a Neumann condition like

$$ D(u^{(n)}) \frac{\partial u^{(n+1)}}{\partial x} = F(u^{(n)})$$

where $u^{(n)}$ is the approximation of $u$ at time $t$. You can also try the following Robin condition:

$$ D(u^{(n)}) \frac{\partial u^{(n+1)}}{\partial x} = u^{(n+1)} v (1-u^{(n)})^{C}$$

This should give you an approximative conservation of the mass (the conservation should also improve as the time step decreases). If you want to conserve exactly the mass, then you should subiterate at each time step to solve the non-linearity (but then the convergence of the non-linear solver comes into play).

Hope this helps!


I guess so. At least i get an output. But somehow my mass is not always conserved (depends on proper initial values, which is bad for data fitting)

I tried today following: I replaced the flux function $F(u) = u*v*(1-u)^C $with $ F(u) = u*v*sign(1-u)*|1-u|^C$

so i can set a boundary like $v*sign(1-u)*|1-u(end)|^C = -v*sign(1-u)*|1-u(end-1)|^C $

so the solver sets the concentration in the last cell > 1 so that everything flowing out is flowing back. That conserves mass but takes longer and sometimes does not converge.

my boundaries right now are

if hyperbolic





Oh, man. Just came over this old thread. If anybody reads this, don´t use toolboxes. Write it your own. I did it with the method of lines and published boundary conditions. Working well with standard ode solvers. Don´t know why it did not work.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.