-- 
_________________________________________________

Timo Koch                                      phone: +49 711 685 64676
IWS, Universität Stuttgart                  fax:   +49 711 685 60430
Pfaffenwaldring 61         email: [email protected]
D-70569 Stuttgart             url: www.iws.uni-stuttgart.de/en/lh2/
_________________________________________________

> On 26. May 2020, at 17:13, Dmitry Pavlov <[email protected]> wrote:
> 
> 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?

Puh on first look I would say you are correct. As the element solution is 
passed to the spatial params interfaces, you would expect that you are allowed 
to use all of the values without messing with the Jacobian. So this would 
clearly be a bug in Dumux / in the numeric differentiation assembler.
It is probably assumed that the volume variables will only depend on the DOF of 
the sub-control volume / adjacent vertex. Updating all volume variables would 
mean a significant overhead, though.
I’ll have to look at this more carefully to confirm.

> 
> 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?

Not really. You’d have to modify the local assembler to also update all other 
volume variables in the element whenever one of the solutions is deflected. You 
could hack it in there and see if it works as expected.

> 
> 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?

If above is true, a hand-made analytic Jacobian should not be affected by this 
problem, only the assembler using numeric differentiation.


I hope, we can clear this up on our next core developer meeting (which should 
be tomorrow if I’m not mistaken).

Best wishes
Timo


> 
> (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

_______________________________________________
Dumux mailing list
[email protected]
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux

Reply via email to