How do we correctly define the flux in a finite volume method applied to Poisson's equation where we have a piecewise constant material? Specifically, say we have the equation \begin{align*} -\nabla\cdot(k\nabla u) =& f\textrm{ on }\Omega\\ u =& 0\textrm{ on }\partial\Omega \end{align*} where $\Omega$ is some kind of easy, square domain in $\mathbb{R}^2$. Then, let's mesh $\Omega$ using square elements where one of the cells is $K$. Then, on the interior, we can integrate over $K$, which gives $$ \int_{K}-\nabla\cdot(k\nabla u) = \int_{K}f $$ After using the divergence theorem, we have $$ \int_{\partial K}-k(\nabla u\cdot n) = \int_{K}f $$ Now, we have to discretize the quantity on the left, which can be divided into four separate integrals. Assume that center of cell $K$ is $(x_i,y_j)$ and let $u_{ij}=u(x_i,y_j)$ and $k_{ij}=k(x_i,y_j)$. Then, we can discretize $\nabla u \cdot n$ into four different pieces, depending on what side of the square we're on, as \begin{align*} &\frac{u_{ij}-u_{{i-1},j}}{\Delta x}\\ &\frac{u_{ij}-u_{{i+1},j}}{\Delta x}\\ &\frac{u_{ij}-u_{i,{j-1}}}{\Delta y}\\ &\frac{u_{ij}-u_{i,{j+1}}}{\Delta y} \end{align*} What I'm not clear about is what to do about $k$ if we assume that $k$ is piecewise constant and defined on the same mesh as $u$. We could just multiply the fluxes by $k_{ij}$. For example, discretize $-k(\nabla u\cdot n)$ along the left face as $$ k_{ij}\left(\frac{u_{ij}-u_{{i-1},j}}{\Delta x}\right) $$ However, in this case, the flux on the other side would be $$ k_{{i-1},j}\left(\frac{u_{{i-1},j}-u_{i,j}}{\Delta x}\right) $$ so they wouldn't match and be conservative because they're using two different material values. The best I can think of is to average the material and use $$ \left(\frac{k_{ij}+k_{{i-1},j}}{2}\right)\left(\frac{u_{ij}-u_{{i-1},j}}{\Delta x}\right) $$ but I'm not sure if this is correct or the implications of this choice.

As such, how do we discretize the flux when we have a piecewise constant material, $k$, defined on the same mesh as our material $u$.


Taking the average will certainly work, but is recommended to take the harmonic mean to account for the changing conductivity:

$k = \Big( \frac{1-f}{k_i} + \frac{f}{k_{i+1}} \Big)^{-1}$

If you choose your faces directly midway between your two cells, f = 0.5:

$ k = \frac{2k_ik_{i+1}}{k_i + k_{i+1}} $

It can be derived by looking at the physical flux, which you want to approximate.

I suggest, you take a look at the book by Patankar: Numercial Heat Transfer and Fluid Flow Chapter 3. It can be found in any decent science library and I guess, the internet counts as well.


Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

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