Dear Timo, Dear Leo,

Thank you very much for your replies, which are very helpful for me to diagnose 
my DuMuX implementation. I have a few follow-up questions as well. Thank you 
again for your time.


Hi Timo,

Thank you for explaining the default residual/Jacobian evaluation workflow. 
Yes, my model only couples via source terms, and I use the BOX discretization. 
Based on your explanation, I believe my current issues are likely not caused by 
the two Jacobian-related questions themselves. I really appreciate the 
clarification, which helped me narrow down the problem significantly.

Regarding the point about the source term being “volume-specific”: in my 
implementation, I compute q_mass = rho * alpha_p * (p_other - p_this), where 
alpha_p * (p_other - p_this) has units of [1/s], so q_mass has units of 
[kg/(m^3·s)]. Under this interpretation, the source I return is already a mass 
source density (per volume), and after local integration it is multiplied by 
the SCV volume to become mass-based, which I believe is consistent with the 
“volume-specific” convention. Therefore, I assumed I would not need to divide 
by scv.volume() again. If I misunderstood what you meant by “volume-specific” 
(e.g., if Problem::source(...) is not expecting a volumetric source density 
with units of [kg/(m^3·s)] in this context), I would be very grateful if you 
could point me to the correct convention for Problem::source(...).

Finally, a question about the “single-domain extension” approach you mentioned: 
for a fully overlapping dual-continuum model with the same discretization, I 
understand this as implementing both continua within one domain by adding an 
additional pressure equation, and including the exchange term as an internal 
coupling between the equations. If that is the right interpretation, could you 
briefly indicate what additional development work this typically requires 
compared to using the multidomain framework? In particular, I am not sure 
whether extending the Richards model to two coupled equations in DuMuX is 
straightforward, or which components would usually need to be modified.


Hi Leo,

Thank you for pointing out the distinction between dual-porosity and 
dual-permeability. Sorry for the ambiguity in my earlier wording: both the 
DuMuX model I built and the HYDRUS-1D model I compared against are 
dual-permeability models (flow is considered in both domains), not 
dual-porosity.

You mentioned your dual-permeability implementation work from 2010. I am not 
sure whether that work or the code was formally published, but I would be very 
interested in reading any related code or materials you are comfortable 
sharing, and I would be happy to cite them appropriately in my ongoing work. If 
any version of the code still exists and can be shared for reference, it would 
be extremely helpful for me to cross-check conventions and implementation 
details while I diagnose my current issue.


I sincerely appreciate both of your time and guidance. If you would like me to 
provide additional details, please let me know. It is exciting to discuss and 
learn DuMuX with both of you. Thank you again for maintaining the project and 
for your kind support.

Best regards,
Wencheng Qi, PhD
Postdoc Research Scientist
Department of Crop Sciences
University of Illinois Urbana–Champaign


From: DuMux <[email protected]> on behalf of 
[email protected] <[email protected]>
Date: Tuesday, January 13, 2026 at 5:30 AM
To: DuMux User Mailing List <[email protected]>, Koch, Timo 
<[email protected]>
Subject: Re: [DuMux] Questions on developing a dual-porosity model in DuMuX 
MultiDomain

Dear Wencheng,

Timo has already mentioned all the important details. So, I can only add a few 
minor points.

If flow occurs in the soil and the macropore domain, the model is called a 
dual-permeability model. In a dual-porosity model, only flow in the macropores 
is considered. Perhaps this explains the observed differences between your 
DuMux dual-permeability model and the HYDRUS-1D dual-porosity model.

I had implemented a two-phase (water/air) dual-permeability model in DuMux in 
2010 (this model does no longer exist in DuMux). Therefore I extended the 
equations and defined the soil parameters for each domain (as Timo mentioned). 
The domains where coupled via the source term, which worked very well. I also 
conducted a comparison between a 1D dual-permeability model based on the 
Richards equation and the two-phase water/air counterpart 
(dual-podual-pohttps://doi.org/10.1007/s11242-012-0064-3<https://urldefense.com/v3/__http://dx.doi.org/10.14279/depositonce-4972__;!!DZ3fjg!7EwLRt5qLFd5w4KcfvhcRIkuivAmwwYXnFGMDPhqB94EHLGTNGKts-ZwYZCYAGT69j7oefQwD_OT9lc7dqSDjlk3FZ4RwQ$>
 , 
http://dx.doi.org/10.14279/depositonce-4972<https://urldefense.com/v3/__http://dx.doi.org/10.14279/depositonce-4972__;!!DZ3fjg!7EwLRt5qLFd5w4KcfvhcRIkuivAmwwYXnFGMDPhqB94EHLGTNGKts-ZwYZCYAGT69j7oefQwD_OT9lc7dqSDjlk3FZ4RwQ$>).
 Perhaps the results of these studies could be used for validation or 
comparison purposes.

Best regards
Leo
Koch, Timo <[email protected]> hat am 13.01.2026 10:38 CET 
geschrieben:


Dear Wencheng,

welcome to the mailing list! I’m glad you found DuMux for this task and find it 
useful.
I’ll try to answer your questions below inline.

On 8 Jan 2026, at 21:03, Qi, Wencheng <[email protected]> wrote:

Dear DuMuX Developer Team,

Happy New Year! Thank you for your continued maintenance and development of 
DuMuX. My name is Wencheng Qi, and I am currently a postdoctoral researcher at 
the University of Illinois Urbana–Champaign working on soil–plant interactions. 
I am developing a dual-porosity model using the DuMuX Multidomain framework, 
and I have some questions regarding Jacobian assembly and cross-domain 
derivatives in Multidomain. Thank you very much for your time, and I would 
greatly appreciate any guidance you could provide.

The dual-porosity model I am trying to implement is conceptually similar to the 
root–soil two-domain coupling approach in dumux/test: the two continuum domains 
represent the matrix and macropore regions, and each domain solves a Richards 
equation with its own Problem. Unlike an explicit geometric partitioning (e.g., 
“fracture” vs “matrix”), my two domains are conceptually overlapping: at each 
spatial location, both matrix and macropore states exist, and mass is exchanged 
between the two domains via an exchange term.

If the the domain are fully overlapping, so essentially you have a true 
dual/multi continuum model, and your are also using the same discretization 
method for all equations, it is usually more efficient and effective to extend 
a given model (e.g. Richards) by more equations rather than realizing this in a 
multi-domain framework. The reason is that the Jacobian derivatives can be 
computed a bit more efficiently in this case. Also the solver can utilize the 
block structure. However, this is not to say that it doesn’t work with 
multi-domain. With multidomain you get the Jacobian in global model-specific 
block structure which might also be useful to you.


After the model development, I tried to conduct a model validation between 
DuMuX and HYDRUS-1 (with dual-porosity). However, there is no agreement in flow 
variables (e.g., water exchange flux) between the two models. Thus, I have two 
major questions:


  1.
After adding the exchange term as a source term to the governing equation in 
each Problem, will this contribution automatically be included in 
evalLocalResidual? In particular, when DuMuX builds the Jacobian by numerical 
differentiation, will the newly added source term be differentiated and 
assembled into the Jacobian automatically? In other words, in typical DuMuX 
workflows, do I not need to manually implement the Jacobian specifically for 
this added source term


Yes this is automatic (using numeric differences). In case you haven’t 
overloaded the evalCouplingResidual function, per default, the full residual 
will be evaluated which in includes all terms.


  1.
Because the exchange term depends on both the pressure in the current domain 
(p_this) and the pressure in the other domain (p_other), the Jacobian should, 
in principle, contain cross-domain block entries. Will these cross-domain 
derivatives also be included automatically? Or is it recommended to implement 
the exchange via a specific coupling residual/interface (rather than directly 
adding it in each Problem’s source(...)), or to explicitly provide/supplement 
Jacobian contributions somewhere in the CouplingManager/Assembler?



Yes, the cross-domains are also computed automatically based on what you 
provide in “couplingStencil”. The couplingStencil has to return a list of all 
degrees of freedom of the other domain that have an influence on the local 
residual of the given element in this domain. Often, if you forget some degrees 
od freedom, you will obtain only a rough approximation to the Jacobian and this 
increases the number of Newton iteration but you can often still get 
convergence. However, this is not guaranteed. Note that the dofs you return 
depend on the discretization scheme. For BOX/CVFE in your case it would 
typically be all dofs of the element, whereas for CC-TPFA/MPFA it would be the 
element dof. In case your coupling stencil only contains the dof of the same 
control volume for TPFA, only the source/storage derivatives will be evaluated 
correctly, not the fluxes, but I assume your model only couples via source 
terms.



Supplementary details

My implementation workflow is as follows: I created two Richards subproblems 
(problem-matrix and problem-macropore) and implemented a CouplingManager to 
handle the mapping and exchange computation between the two domains. In the 
CouplingManager, I define a one-to-one DoF mapping (each grid cell corresponds 
to the cell at the same location in the other domain), and provide an 
exchangeSource(domain, scv) interface. In that function, I use the mapping to 
find the DoF in the other domain, read the pressure variables from curSol, 
compute the exchange term in the form q_mass = rho * alpha_p * (p_other − 
p_this) (where alpha_p is computed from my spatial parameters/hydraulic 
conductivity/geometric factors, etc.), and return the mass source density 
q_mass. Then, in each Problem I implement source(...) and add 
couplingManager().exchangeSource(...) to the continuity equation.

My current understanding is: once the exchange term is added via source(...), 
it becomes one of the contributions to the local residual together with the 
storage term and the flux term. Therefore, when evalLocalResidual evaluates the 
total residual, it should evaluate these contributions together (storage + flux 
+ source). If the Jacobian is then built by numerical differentiation, 
DuMuXperturbs the solution vector and re-evaluates the same total residual; in 
that process, the exchange source is re-evaluated as well and thus its 
derivatives should be captured automatically in the numerically differentiated 
Jacobian. Under this understanding, I would not need to add a separate Jacobian 
implementation specifically for the exchange source term.

Yes this is correct.


A simplified version of my CouplingManager exchangeSource is:

// CouplingManager: compute exchange using mapped dofs and current solutions
Scalar exchangeSource(domainI, scv) const
{
dofThis = scv.dofIndex();
dofOther = mapThisToOther[dofThis];
pThis  = curSol(domainI)[dofThis ][pIdxThis ];
pOther = curSol(other )[dofOther][pIdxOther];
return q_mass = rho * alpha_p * (pOther - pThis);
}


And the source implementation in problem-macropore is:

PrimaryVariables source(const Element& element,
const FVElementGeometry& fvGeometry,
const ElementVolumeVariables& elemVolVars,
const SubControlVolume& scv) const
{
PrimaryVariables src(0.0);
const Scalar q_mass = 
couplingManager().exchangeSource(Dune::index_constant<1>{}, scv);
src[Indices::conti0EqIdx] += q_mass;
return src;
}


The source in DuMux should be volume-specific. This might be the issue here. If 
you have the full mass here, you need to divide by the scv.volume() to have a 
volume-specific term.
(The term gets integrated and thus multiplied by the volume again later in the 
local residual. But commonly/historically volume-specific source have been a 
bit more convenient, not here.)

Best wishes and let us know if you have any follow-up questions.
Timo

P.S. You can also reach out to us through the matrix channel.

——
__________________________________________________________

Jun.-Prof. Dr.-Ing. Timo Koch

Department for Hydromechanics and Modelling of Hydrosystems (LH2)
Institute for Modeling Hydraulic and Environmental Systems (IWS)
University of Stuttgart, Germany

Institute webpage:   
https://www.iws.uni-stuttgart.de/en/lh2/<https://urldefense.com/v3/__https://www.iws.uni-stuttgart.de/en/lh2/__;!!DZ3fjg!7EwLRt5qLFd5w4KcfvhcRIkuivAmwwYXnFGMDPhqB94EHLGTNGKts-ZwYZCYAGT69j7oefQwD_OT9lc7dqSDjlndBjd_uA$>
Personal webpage: 
https://timokoch.github.io/<https://urldefense.com/v3/__https://timokoch.github.io/__;!!DZ3fjg!7EwLRt5qLFd5w4KcfvhcRIkuivAmwwYXnFGMDPhqB94EHLGTNGKts-ZwYZCYAGT69j7oefQwD_OT9lc7dqSDjln2tMSRIw$>
Phone:                    +49 711 685 60461
___________________________________________________________



If anything in my description or understanding is incorrect, I would greatly 
appreciate your corrections. I am also happy to provide additional 
implementation details to help you diagnose the issue.

Thank you very much!

Best regards,
Wencheng Qi, PhD
Postdoc Research Scientist
Department of Crop Sciences
University of Illinois Urbana–Champaign
_______________________________________________
DuMux mailing list
[email protected]
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux<https://urldefense.com/v3/__https://listserv.uni-stuttgart.de/mailman/listinfo/dumux__;!!DZ3fjg!7EwLRt5qLFd5w4KcfvhcRIkuivAmwwYXnFGMDPhqB94EHLGTNGKts-ZwYZCYAGT69j7oefQwD_OT9lc7dqSDjlnTn2XuSA$>
_______________________________________________
DuMux mailing list
[email protected]
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux<https://urldefense.com/v3/__https://listserv.uni-stuttgart.de/mailman/listinfo/dumux__;!!DZ3fjg!7EwLRt5qLFd5w4KcfvhcRIkuivAmwwYXnFGMDPhqB94EHLGTNGKts-ZwYZCYAGT69j7oefQwD_OT9lc7dqSDjlnTn2XuSA$>


Im Auftrag

Dr.-Ing. Leopold Stadler

--
Referat Daten und Methoden
Abteilung Wasserbau im Binnenbereich

Bundesanstalt für Wasserbau
Federal Waterways Engineering and Research Institute
Kußmaulstraße 17 | 76187 Karlsruhe

https://www.baw.de<https://urldefense.com/v3/__https://www.baw.de__;!!DZ3fjg!7EwLRt5qLFd5w4KcfvhcRIkuivAmwwYXnFGMDPhqB94EHLGTNGKts-ZwYZCYAGT69j7oefQwD_OT9lc7dqSDjlkmtVoctQ$>

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

Reply via email to