Hi Johannes, Thank you for your comments.
Looking into the issue some more detail, it seems our fluxes need to be calculated in mass terms because the FluxVariables property has been set to the ImplicitDarcyFluxVariables template instead of the TwoPNCFluxVariables template which uses mole fraction gradients and the DiffusivityMillingtonQuirk model. I have changed our boundary and initial conditions to molar terms, converting the phase flow which contains the minor components to mole flow in the problem template using the average molecular weight, and the results look good to me at the moment. Meanwhile I am consulting with Dr. Martin Diaz who is more familiar with the results which should be expected. Kind regards, Edscott -----Mensaje original----- De: Johannes Hommel [mailto:[email protected]] Enviado el: jueves, 28 de septiembre de 2017 03:51 a. m. Para: Ed Scott Wilson Garcia; DuMuX User Forum Asunto: Re: [DuMuX] error in 2pnc model? Hello, you are probably right that the difference results from the different formulations (mass/mole). For the compositional models, we use mole fractions as primary variables as this facilitates the calculation of correct diffusive fluxes, phase transitions by fugacities and reactive processes. And while it is true that the overall mass is conserved during reactions and the number of moles isn't necessarily, the neither mass or mole fraction of a components are conserved when it is involved in reactions. To get the same results as for the 2p case, you would need to convert your 2p case boundary conditions from mass to mole-based quantities, so that the conditions are in the same units as the flux calculations in the 2pnc model. Can you try that and see whether it fixes your problem? Best regards, Johannes On 09/22/2017 08:38 PM, Ed Scott Wilson Garcia wrote: > 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 -- ******************************************************* Johannes Hommel Lehrstuhl für Hydromechanik und Hydrosystemmodellierung (LH2) Institut für Wasser- und Umweltsystemmodellierung (IWS), Universität Stuttgart www.hydrosys.uni-stuttgart.de email: [email protected] Pfaffenwaldring 61 70569 Stuttgart Tel.: ++49 711 / 685-64600 Department of Hydromechanics and Modelling of Hydrosystems (LH2) Institute for Modelling Hydraulic and Environmental Systems (IWS), University of Stuttgart www.hydrosys.uni-stuttgart.de email: [email protected] Pfaffenwaldring 61 70569 Stuttgart Tel.: ++49 711 / 685-64600 _______________________________________________ Dumux mailing list [email protected] https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
