GradientFluxAdapter problem
Incorrect representation of fluxes in the output
Summary
Since fluxes are quite important for the transport model, in dune-modelling
they are calculated directly from the local operator, later reconstructed and stored using Raviart-Thomas elements, meanwhile, in DORiE
they are calculated only for postprocessing purposes using GridFunctionInterface
and later adapted to the vtk output with VTKGridFunctionAdapter
. Details of both are quite diferent, but in the end, I would say they must representing the very same flux!
Steps to reproduce
Enforce a strong Neumann boundary condition.
What is the current bug behaviour?
Here I'm posting the differences I have found between both for a 1D case with infiltration. (left DORiE
VTKAdapter, right dune-modelling
flux reconstruction)
In DORiE
, the flux is being calculated point-wise, which is wrong because in the DG model they have only meaning in the interface of elements. Moreover, the evaluation of the fluxes in the interface is also wrong because they are being computed without taking into account that the head is a discontinuous quantity. This can be seen in the figure: Although the influx is 5.5e-8, only the "inner" flux of the element is computed in the DORiE
adapter. The missing part of the flux is given by the jump between the boundary and the element face. In this case, with a big jump due to a strong boundary condition, the output already differs by three orders of magnitude.
What is the expected correct behaviour?
If the printing is continuous, then, fluxes must be reconstructed. Another option would be to enforce printing only vertex values in the interfaces and restrict somehow the interpolation in paraview such that is clear it doesn't have a clear continuous representation.