# The effect of grid size on the total flux when solving Darcy flow with mixed finite element method

I am solving Darcy flow now with mixed finite element method. The Dary flow is \begin{aligned}k^{-1}\mathbf{q} + \nabla h=0, \text{ in } \Omega\\ % \nabla\cdot \mathbf{q} = 0, \text{ in } \Omega, \\ % h =h_D, \text{ on } \partial\Omega_D, \\% \mathbf{q}\cdot \mathbf{n} = q_N = 0, \text{ on } \partial\Omega_N,% \end{aligned} where $$\mathbf{q}$$ is the flux velocity, and $$h$$ is hydraulic head, $$\mathbf{n}$$ is the normal direction of the boundary.

The boundary conditions are like the figure below.

The mixed FEM method is just a FEM that determine the normal flux velocity on edge and mean head of an element, like the figure below. And I use $$RT0$$ element.

Then, I just calculate the total flux through the domain. To be more exact, the flux enters at the top boundary, and leave at the bottom boundary. The totoal flux is calculated by summing up the product of the flux velovity and associated edge length, which reads: $$$$Q_{in} = \sum\limits_{i = 1}^{NB}q_{et, i} \times L_{et, i},$$,$$ where $$NB$$ is the number of element edges adhering to the top boundary, $$q_{et, i}$$ is the normal flux velocity at the edge, and $$L_{et, i}$$ is the edge length.

If the model domain is regular, like a square/rectange, the grid size does not affect the total flux value.

But, if the domain is not regular, e.g. trapzoidal domain, like the figrue below.

The total flux is significantly influenced by the grid/cell size. The relation is, the grid size is smaller, the total flux is higher, like the plot below.

My questions are

1. Why the totoal flux increases with the decrease in grid size? why this happens?

2. If the grid size is nearly zero, is that mean: the total flux will be the true value?

## 1 Answer

All you are guaranteed by mathematical analysis is that the numerical solution converges as the mesh size decreases. In some cases, quantities will converge from below, in others from above, and in yet others, convergence may be oscillatory. What it is in your case is not obvious and your observation that convergence appears to be from below seems as good a conjecture as one can come up with short of actually being able to prove it.

Convergence also means that if you choose very small mesh sizes, then the flux value you compute is close to the exact value. So the answer to your question 2 is "yes".