Working on a two liquid immiscible phase, n component problem we seem to have
stumbled across what seems like an error within the 2pnc model. I explain.
We started from the GeneralizedDirichletProblem of the 2p model. After
switching to the pnsw formulation, defining our grid and setting the
appropriate boundary and initial conditions, we obtained satisfactory results.
We followed up by replacing the 2p model for the 2pnc and making the two
fluids (water, lnapl) immiscible with Johannes' suggestion (provided on this
mail list). Setting the third component concentration to zero should give
similar
results to the 2p case. This does not happen.
The transport equation for each component is expressed in density and mass
fraction terms, which is a mass tranfer equation. Nonetheless, in
2pnc/implicit/localresidual.hh the advective flux and storage terms are
calculated using the molar density and molar fraction, which yields a quantity
representing a number of atoms divided by Avogadro's number. Although mass is
always conserved, this is not necessarily the case with moles, as the number
of molecules is subject to the chemical reactions which may occur.
Converting the advective and storage terms to mass is easy, just multiply by
the component's molecular mass. By doing so (and adding the necessary upstream
methods to the compositional fluid state) the results obtained are comparable
to the 2p case. Below you will find diff describing the modifications (extra
parenthesis are not necessary, but I like them).
Is there any reason why the advective and storage terms in the 2pnc model are
being calculated in number of atoms instead of mass? Would this be some
particularity of the fluidcellproblem?
Any comments will be greatly appreciated.
kind regards,
Dr.Edscott Wilson Garcia
Reservoir Engineering
Mexican Petroleum Institute
--- GIT/dune/dumux/dumux/porousmediumflow/2pnc/implicit/localresidual.hh
2017-03-31 12:40:09.000000000 -0600
+++ GIT/LSWF/tests/include/2p3c-3/localresidual.hh 2017-09-21
13:06:10.000000000 -0500
@@ -176,16 +176,18 @@
int eqIdx = conti0EqIdx + compIdx;
if (replaceCompEqIdx != eqIdx)
{
- storage[eqIdx] += volVars.molarDensity(phaseIdx)
+ storage[eqIdx] += (volVars.molarDensity(phaseIdx)
* volVars.saturation(phaseIdx)
* volVars.moleFraction(phaseIdx, compIdx)
- * volVars.porosity();
+ * volVars.porosity()
+ * volVars.molarMass(compIdx));
}
else
{
- storage[replaceCompEqIdx] += volVars.molarDensity(phaseIdx)
+ storage[replaceCompEqIdx] +=
(volVars.molarDensity(phaseIdx)
* volVars.saturation(phaseIdx)
- * volVars.porosity();
+ * volVars.porosity()
+ * volVars.molarMass(compIdx));
}
}
}
@@ -242,17 +244,19 @@
if (massUpwindWeight_ > 0.0)
// upstream vertex
flux[eqIdx] +=
- fluxVars.volumeFlux(phaseIdx)
+ (fluxVars.volumeFlux(phaseIdx)
* massUpwindWeight_
* up.molarDensity(phaseIdx)
- * up.moleFraction(phaseIdx, compIdx);
+ * up.moleFraction(phaseIdx, compIdx)
+ * up.molarMass(compIdx));
if (massUpwindWeight_ < 1.0)
// downstream vertex
flux[eqIdx] +=
- fluxVars.volumeFlux(phaseIdx)
+ (fluxVars.volumeFlux(phaseIdx)
* (1 - massUpwindWeight_)
* dn.molarDensity(phaseIdx)
- * dn.moleFraction(phaseIdx, compIdx);
+ * dn.moleFraction(phaseIdx, compIdx)
+ * dn.molarMass(compIdx));
Valgrind::CheckDefined(fluxVars.volumeFlux(phaseIdx));
Valgrind::CheckDefined(up.molarDensity(phaseIdx));
@@ -265,15 +269,17 @@
// upstream vertex
if (massUpwindWeight_ > 0.0)
flux[replaceCompEqIdx] +=
- fluxVars.volumeFlux(phaseIdx)
+ (fluxVars.volumeFlux(phaseIdx)
* massUpwindWeight_
- * up.molarDensity(phaseIdx);
+ * up.molarDensity(phaseIdx)
+ * up.molarMass(compIdx));
// downstream vertex
if (massUpwindWeight_ < 1.0)
flux[replaceCompEqIdx] +=
- fluxVars.volumeFlux(phaseIdx)
+ (fluxVars.volumeFlux(phaseIdx)
* (1 - massUpwindWeight_)
- * dn.molarDensity(phaseIdx);
+ * dn.molarDensity(phaseIdx)
+ * dn.molarMass(compIdx));
Valgrind::CheckDefined(fluxVars.volumeFlux(phaseIdx));
Valgrind::CheckDefined(up.molarDensity(phaseIdx));
Valgrind::CheckDefined(dn.molarDensity(phaseIdx));
@@ -328,10 +335,11 @@
result = 0;
for (int compIdx = 0; compIdx < numComponents; ++compIdx)
{
- result[conti0EqIdx + compIdx] += volVars.density(phaseIdx)
+ result[conti0EqIdx + compIdx] += (volVars.density(phaseIdx)
*
volVars.saturation(phaseIdx)
*
volVars.massFraction(phaseIdx, compIdx)
- * volVars.porosity();
+ * volVars.porosity()
+ * volVars.molarMass(compIdx));
}
result *= this->fvGeometry_().subContVol[scvIdx].volume;
}
_______________________________________________
Dumux mailing list
[email protected]
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
