Hello,
I am trying to do a 1D porous medium flow simulation where one of the
components is a surfactant that affects krw and krn. The effect of the
surfactant depend on pressure gradient.
Earlier, an attempt [1] was made to apply CC method for this kind of
problem. It turned out that DuMux API does not easily allow to estimate
a gradient in TPFA, and, following Timo Koch's advice, I gave a try to
Box method.
I think I am having some trouble with derivatives now. Looking at this
comment in boxlocalassembler.hh, I am beginning to understand why.
// Calculate derivatives of all dofs in stencil with respect to
the dofs in the element. In the //
// neighboring elements we do so by computing the derivatives
of the fluxes which depend on the //
// actual element. In the actual element we evaluate the
derivative of the entire residual. //
Why the trouble? Well, I am calculating the pressure gradient (by
calling evalGradients) inside this method
MaterialLawParams materialLawParams(const Element& element,
const SubControlVolume& scv,
const ElementSolution& elemSol)
const
Here, I store the needed numbers for krw and krn calculation in
MaterialLawParams, and they are properly downstream at MaterialLaw::krw
and MaterialLaw::krn.
Now, let my have two neighbor boxes, 0 and 1. Pressure in box 1 affects
the krw/krn in box 1. evalGradients, when called, duly calculates
pressure gradient in box 0 depending, among others, on the pressure
value in box 1. That is good.
But it turns out that the numerical differentiation algorithm does not
bother to call materialLawParams() for scv-s in box 0 when it calculates
the derivatives of fluxes in box 0 w.r.t. pressure in box 1. I suppose
that it has to do with the comment above. I suppose that it takes into
account only the "transmissibility" part of the effect of pressure in
box 1 to flux in box 0, and skips the part that comes from krw/krn
sensitivity to the pressure gradient.
Also this comment in boxlocalassembler.hh may be relevant too.
// TODO additional dof dependencies
I will very much appreciate answers to the following questions:
1. Am I correct about the behavior of the numerical Jacobian assembler,
or it should be all right and there is a bug somewhere in my code or
DuMux's?
2. In case I am correct, is there an easy (or not) way to force DuMux
into recalculating the material law parameters for scv-s that belong to
DOFs w.r.t. which we take the derivative? Or take some other approach?
3. If not, will the hand-made analytic Jacobian help? Or it will fail
due to some other assumption in the DuMux engine that is not true for my
problem?
(I actually have the analytic Jacobian implemented, but not quite sure
it is correct, partially because, well, I never had a chance to test it
against a correct numeric Jacobian).
Thank you for your time.
Best regards,
Dmitry
[1] https://listserv.uni-stuttgart.de/pipermail/dumux/2020q2/002516.html
_______________________________________________
Dumux mailing list
[email protected]
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
