# Discretizing a parabolic PDE with finite volume method

I want to discretize the following parabolic PDE:

$$u_t = \nabla\cdot(\alpha(x)\nabla u)- \beta u\\ x\in\Omega \subset \mathbb{R}^2\\ \partial_n u = 0\\ u(t,0) = u_0(x)\ge 0, \alpha(x)>0$$

Given the Neumann boundary condition above, one can integrate this PDE:

\begin{alignat}{2} \frac{\mathrm d}{\mathrm dt}\int_\Omega \! u\,\mathrm dx&=\phantom{-}\int_\Omega \! \nabla\cdot (a(x)\nabla u)\,\mathrm dx \:\:&-\int_\Omega \! \beta(x) u \,\mathrm dx \\ &=\phantom{-}\int_{\partial\Omega}\! a(x)\nabla u\cdot n \,\mathrm ds&-\int_\Omega \! \beta(x) u \, \mathrm dx \\ &=\phantom{-}\int_{\partial\Omega} a(x)\partial_n u \,\mathrm ds&-\int_\Omega \! \beta(x) u \, \mathrm dx \\ &=-\int_\Omega \! \beta(x) u \, \mathrm dx \end{alignat}

The form of the discretized linear system I want to get is

$$\frac{d\vec{u}}{dt}=A\vec{u}-\vec{c}$$

Unfortunately, I don't think this is exactly what I'm obtaining, despite I do think I use a correct discretization technique. To this end, please consider my steps below.

I discretize the two integrals using the finite volume method:

$$\frac{d}{dt}\sum\limits_{j\in n_i} u_j l_{ij}=-\sum\limits_{j\in n_i} \beta_{ij} u_j l_{ij},$$

where $$n_i$$ is the set of neighboring cells to the cell $$V_j\subset \Omega$$ and where $$l_{ij}$$ is the length of the boundary of $$V_j$$ between $$u_i$$ and $$u_j$$.

Now one can discretize it in time to obtain:

$$\sum\limits_{j\in n_i} \frac{u_j^{n+1}-u_j^n}{\Delta t}l_{ij}=-\sum\limits_{j\in n_i} \beta_{ij} u_j^n l_{ij}$$

So we obtain the linear system $$L\frac{d\vec{u}}{dt} = -A\vec{u}+\vec{c},$$

where $$L$$ is the matrix containing $$l_{ij}$$ and each entry of $$A$$ contains $$-\beta_{ij}l_{ij}$$, and where $$\vec{c}$$ contains terms of the boundary condition.

I would appreciate if someone could please let me know if this is correct. I think I'm quite confused as I don't see where each term of the matrix should go, and where the vector $$\vec{c}$$ is supposed to be.

Does this mean that the correct equation in the form $$\frac{du}{dt}=Ab+c$$ would be, in this case,

$$\frac{d\vec{u}}{dt} = -L^{-1}A\vec{u}+L^{-1}\vec{c}$$

• Think over which domain you integrate. You should integrate over the subdomain of the entire domain, a d thus you cannot apply the boundary conditions given at the boundary of the global domain – VorKir Mar 28 at 3:29
• @VorKir Thanks but I don't have that much insight to even apply your hint. – sequence Mar 29 at 8:37
• You can find details in many places, e.g., the first link in my Google search google.com/url?sa=t&source=web&rct=j&url=http://… – VorKir Mar 30 at 1:33