I am trying to solve this diffusion equation : $\dfrac{\partial D\dfrac{\partial f}{\partial x}}{\partial x}+S = \dfrac{\partial f}{\partial t}$ ($D$ is not constant and varies according to $x$) with the following BC: $f(x,0)=1 , f(0,t)=0, \dfrac{\partial f}{\partial x}(1,t)=0$

I am using a central finite difference scheme (2nd order) for space and Euler explicit for time (1st order). The discretized $[0..1]$ domain contains the $x_{i}$ points, $i \in [0..N]$ ($x_{0}=0 , x_{N}=1$).

The implementation of Dirichlet condition at $x=0$ is simple. For the Neumann boundary condition, I saw a lot of references that treated simple cases, adequate for constant diffusion coefficients and stationnary cases, like the ghost point, that I can't consider in my case(Diffusion coefficient can't be evaluated out of the domain, ghost point envolves evaluation the equation on a point that is not inside the interior domain... ).

So I am trying something at my own and I want confirmation or correction if I omitted an important aspect in my work:

The discretized equation evaluated at $x_{N-1}$ involves the value of $f_{N}$. If we had Diricihlet condition at $x_{N}$, the problem would have been solved directly. But in our case, $f_{N}$ is unknown.

What I did is that I treated the Neumann BC with a backward scheme (at the 2nd order to preserve the general order of the whole scheme) in order to have an expression of $f^{k}_{N}$:

$\dfrac{\partial f}{\partial x}(x_{N},k\Delta t)=\dfrac{-3f^{k}_{N}+4f^{k}_{N-1}-f^{k}_{N-2}}{2\Delta x}=0$

$f^{k}_{N}=4/3 f^{k}_{N-1}-1/3 f^{k}_{N-2}$

then I injected the new expression of $f^{k}_{N}$ in the discretized equation at $x_{N-1}$. Meaning that I did not add an equation to the problem as I would do for a stationnary problem but I injected the expression of $f^{k}_{N}$ inside the last equation.

So the unknown vector $f^{k}$ in a case of matrix writing would only contain $f^{k}_{i}, i\in [1..N-1]$ (like in a full Dirichlet situation)

Is it correct ?

For information, the matrix form of the discretized problem is : $f^{k+1}=Af^{k}+S+B$, with $B$ the vector containing the B.C , and $S$ the 2nd member vector

  • $\begingroup$ Note that I solved the stationary problem adding a line in the resolution matrix related to the unknown $f_{N}$ containing the Backward finite difference 2nd order discretization at the boundary and worked fine. But I couldn't do the same for an unstationary problem. $\endgroup$ Apr 2, 2018 at 12:12
  • $\begingroup$ The treatment of the stationary problem is presented here : scicomp.stackexchange.com/questions/7688/… $\endgroup$ Apr 2, 2018 at 12:52
  • $\begingroup$ Yes It is correct. What you said: "A full dirichlet situation" does not imply that you have a square matrix of size $N-1$. At the end you must have $N+1$ equations, in this case one is imposed at the beginning (Dirichlet) and the other is calculated once the system is solved (Neumann). Therefore it is not a "full Dirichlet situation",otherwise the boundary values would have been known "a priori". $\endgroup$
    – HBR
    Apr 2, 2018 at 15:58
  • $\begingroup$ @HBR Thanks for your reply. What I mean by a full Dirichlet situation is when the Neumann condition is also replaced by a Dirichlet condition (2 Dirichlet condition for the whole problem , at $x=0$ and at $x=1$). In this case $f^{k}_{N} $ is known, and the problem square matrix has a size of $N-1$ $(N+1-2)$. Now in the case of a Neumann boundary for $x_{N}$ we have to add an equation in the last line of the matrix for $f^{k}_{N} $ is unknown and I don't see how to express it in the unstationary case. $\endgroup$ Apr 2, 2018 at 16:13
  • $\begingroup$ This is not a Dirichlet BC is a Neumann BC. You think (wrong) that it is Dirichlet because you can obtain an explicit solution for this node. A Dirichlet BD is an imposed value, which does not depend on other solves variables... $\endgroup$
    – HBR
    Apr 2, 2018 at 16:16


Your Answer

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