Commit 150291e2 authored by Thomas Fetzer's avatar Thomas Fetzer
Browse files

[material,volumevariables]

introduced setMassFraction() in the compositionalfluidstate.hh so that
in case of mass fraction formulation, they can be set directly, the
converting to mole fractions is now handled by the
compositionalfluidstate class

some test failed (kuevette with 5% difference), new references are
committed, however maybe some thresholds have to be adapted

reviewed by kissinger


git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@15433 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent ca3b1347
......@@ -156,18 +156,8 @@ public:
}
else // mass-fraction formulation
{
Scalar massOrMoleFrac[numComponents];
massOrMoleFrac[transportCompIdx] = priVars[massOrMoleFracIdx];
massOrMoleFrac[phaseCompIdx] = 1.0 - massOrMoleFrac[transportCompIdx];
// calculate average molar mass of the gas phase
Scalar M1 = FluidSystem::molarMass(transportCompIdx);
Scalar M2 = FluidSystem::molarMass(phaseCompIdx);
Scalar X2 = massOrMoleFrac[phaseCompIdx];
Scalar avgMolarMass = M1*M2/(M2 + X2*(M1 - M2));
// convert mass to mole fractions and set the fluid state
fluidState.setMoleFraction(phaseIdx, transportCompIdx, massOrMoleFrac[transportCompIdx]*avgMolarMass/M1);
fluidState.setMoleFraction(phaseIdx, phaseCompIdx, massOrMoleFrac[phaseCompIdx]*avgMolarMass/M2);
// setMassFraction() has only to be called 1-numComponents times
fluidState.setMassFraction(phaseIdx, transportCompIdx, priVars[massOrMoleFracIdx]);
}
}
......
......@@ -131,19 +131,16 @@ public:
fluidState.setPressure(phaseIdx, priVars[pressureIdx]);
Scalar x1 = priVars[massOrMoleFracIdx]; //mole or mass fraction of component 1
if(!useMoles) // mass-fraction formulation
if(useMoles)
{
// convert mass to mole fractions
Scalar M0 = FluidSystem::molarMass(phaseCompIdx);
Scalar M1 = FluidSystem::molarMass(transportCompIdx);
//meanMolarMass if x1_ is a massfraction
Scalar meanMolarMass = M0*M1/(M1 + x1*(M0 - M1));
x1 *= meanMolarMass/M1;
fluidState.setMoleFraction(phaseIdx, phaseCompIdx, 1 - priVars[massOrMoleFracIdx]);
fluidState.setMoleFraction(phaseIdx, transportCompIdx, priVars[massOrMoleFracIdx]);
}
else
{
// setMassFraction() has only to be called 1-numComponents times
fluidState.setMassFraction(phaseIdx, transportCompIdx, priVars[massOrMoleFracIdx]);
}
fluidState.setMoleFraction(phaseIdx, phaseCompIdx, 1 - x1);
fluidState.setMoleFraction(phaseIdx, transportCompIdx, x1);
typename FluidSystem::ParameterCache paramCache;
paramCache.updatePhase(fluidState, phaseIdx);
......
......@@ -316,34 +316,17 @@ public:
else if (phasePresence == nPhaseOnly) {
// only the nonwetting phase is present, i.e. nonwetting phase
// composition is stored explicitly.
if(useMoles) // mole-fraction formulation
if(useMoles)
{
Scalar moleFractionN[numComponents];
moleFractionN[wCompIdx] = priVars[switchIdx];
moleFractionN[nCompIdx] = 1 - moleFractionN[wCompIdx];
// set the fluid state
fluidState.setMoleFraction(nPhaseIdx, wCompIdx, moleFractionN[wCompIdx]);
fluidState.setMoleFraction(nPhaseIdx, nCompIdx, moleFractionN[nCompIdx]);
fluidState.setMoleFraction(nPhaseIdx, nCompIdx, 1 - priVars[switchIdx]);
fluidState.setMoleFraction(nPhaseIdx, wCompIdx, priVars[switchIdx]);
}
else // mass-fraction formulation
else
{
// extract _mass_ fractions in the nonwetting phase
Scalar massFractionN[numComponents];
massFractionN[wCompIdx] = priVars[switchIdx];
massFractionN[nCompIdx] = 1 - massFractionN[wCompIdx];
// calculate average molar mass of the nonwetting phase
Scalar M1 = FluidSystem::molarMass(wCompIdx);
Scalar M2 = FluidSystem::molarMass(nCompIdx);
Scalar X2 = massFractionN[nCompIdx];
Scalar avgMolarMass = M1*M2/(M2 + X2*(M1 - M2));
// convert mass to mole fractions and set the fluid state
fluidState.setMoleFraction(nPhaseIdx, wCompIdx, massFractionN[wCompIdx]*avgMolarMass/M1);
fluidState.setMoleFraction(nPhaseIdx, nCompIdx, massFractionN[nCompIdx]*avgMolarMass/M2);
// setMassFraction() has only to be called 1-numComponents times
fluidState.setMassFraction(nPhaseIdx, wCompIdx, priVars[switchIdx]);
}
// calculate the composition of the remaining phases (as
// well as the densities of all phases). This is the job
// of the "ComputeFromReferencePhase" constraint solver
......@@ -400,31 +383,15 @@ public:
// composition is stored explicitly.
if(useMoles) // mole-fraction formulation
{
Scalar moleFractionW[numComponents];
moleFractionW[nCompIdx] = priVars[switchIdx];
moleFractionW[wCompIdx] = 1 - moleFractionW[nCompIdx];
// set the fluid state
fluidState.setMoleFraction(wPhaseIdx, wCompIdx, moleFractionW[wCompIdx]);
fluidState.setMoleFraction(wPhaseIdx, nCompIdx, moleFractionW[nCompIdx]);
fluidState.setMoleFraction(wPhaseIdx, wCompIdx, 1-priVars[switchIdx]);
fluidState.setMoleFraction(wPhaseIdx, nCompIdx, priVars[switchIdx]);
}
else // mass-fraction formulation
{
// extract _mass_ fractions in the nonwetting phase
Scalar massFractionW[numComponents];
massFractionW[nCompIdx] = priVars[switchIdx];
massFractionW[wCompIdx] = 1 - massFractionW[nCompIdx];
// calculate average molar mass of the nonwetting phase
Scalar M1 = FluidSystem::molarMass(wCompIdx);
Scalar M2 = FluidSystem::molarMass(nCompIdx);
Scalar X2 = massFractionW[nCompIdx];
Scalar avgMolarMass = M1*M2/(M2 + X2*(M1 - M2));
// convert mass to mole fractions and set the fluid state
fluidState.setMoleFraction(wPhaseIdx, wCompIdx, massFractionW[wCompIdx]*avgMolarMass/M1);
fluidState.setMoleFraction(wPhaseIdx, nCompIdx, massFractionW[nCompIdx]*avgMolarMass/M2);
// setMassFraction() has only to be called 1-numComponents times
fluidState.setMassFraction(wPhaseIdx, nCompIdx, priVars[switchIdx]);
}
// calculate the composition of the remaining phases (as
// well as the densities of all phases). This is the job
// of the "ComputeFromReferencePhase" constraint solver
......
......@@ -212,14 +212,9 @@ public:
if(useMoles) // mole-fraction formulation
{
// extract _mole_ fractions in the nonwetting phase
Scalar moleFractionN[numComponents];
moleFractionN[wCompIdx] = priVars[switchIdx];
moleFractionN[nCompIdx] = 1 - moleFractionN[wCompIdx];
// set the fluid state
ParentType::fluidState_.setMoleFraction(nPhaseIdx, wCompIdx, moleFractionN[wCompIdx]);
ParentType::fluidState_.setMoleFraction(nPhaseIdx, nCompIdx, moleFractionN[nCompIdx]);
ParentType::fluidState_.setMoleFraction(nPhaseIdx, wCompIdx, priVars[switchIdx]);
ParentType::fluidState_.setMoleFraction(nPhaseIdx, nCompIdx, 1-priVars[switchIdx]);
// TODO give values for non-existing wetting phase
Scalar xwCO2 = FluidSystem::equilibriumMoleFraction(ParentType::fluidState_, paramCache, wPhaseIdx);
......@@ -229,20 +224,8 @@ public:
}
else // mass-fraction formulation
{
// extract _mass_ fractions in the nonwetting phase
Scalar massFractionN[numComponents];
massFractionN[wCompIdx] = priVars[switchIdx];
massFractionN[nCompIdx] = 1 - massFractionN[wCompIdx];
// calculate average molar mass of the nonwetting phase
Scalar M1 = FluidSystem::molarMass(wCompIdx);
Scalar M2 = FluidSystem::molarMass(nCompIdx);
Scalar X2 = massFractionN[nCompIdx];
Scalar avgMolarMass = M1*M2/(M2 + X2*(M1 - M2));
// convert mass to mole fractions and set the fluid state
ParentType::fluidState_.setMoleFraction(nPhaseIdx, wCompIdx, massFractionN[wCompIdx]*avgMolarMass/M1);
ParentType::fluidState_.setMoleFraction(nPhaseIdx, nCompIdx, massFractionN[nCompIdx]*avgMolarMass/M2);
// setMassFraction() has only to be called 1-numComponents times
ParentType::fluidState_.setMassFraction(nPhaseIdx, wCompIdx, priVars[switchIdx]);
// TODO give values for non-existing wetting phase
Scalar xwCO2 = FluidSystem::equilibriumMoleFraction(ParentType::fluidState_, paramCache, wPhaseIdx);
......@@ -273,14 +256,9 @@ public:
if(useMoles) // mole-fraction formulation
{
// extract _mole_ fractions in the nonwetting phase
Scalar moleFractionW[numComponents];
moleFractionW[nCompIdx] = priVars[switchIdx];
moleFractionW[wCompIdx] = 1 - moleFractionW[nCompIdx];
// convert mass to mole fractions and set the fluid state
ParentType::fluidState_.setMoleFraction(wPhaseIdx, wCompIdx, moleFractionW[wCompIdx]);
ParentType::fluidState_.setMoleFraction(wPhaseIdx, nCompIdx, moleFractionW[nCompIdx]);
ParentType::fluidState_.setMoleFraction(wPhaseIdx, wCompIdx, 1-priVars[switchIdx]);
ParentType::fluidState_.setMoleFraction(wPhaseIdx, nCompIdx, priVars[switchIdx]);
// TODO give values for non-existing nonwetting phase
Scalar xnH2O = FluidSystem::equilibriumMoleFraction(ParentType::fluidState_, paramCache, nPhaseIdx);
......@@ -290,20 +268,8 @@ public:
}
else // mass-fraction formulation
{
// extract _mass_ fractions in the nonwetting phase
Scalar massFractionW[numComponents];
massFractionW[nCompIdx] = priVars[switchIdx];
massFractionW[wCompIdx] = 1 - massFractionW[nCompIdx];
// calculate average molar mass of the nonwetting phase
Scalar M1 = FluidSystem::molarMass(wCompIdx);
Scalar M2 = FluidSystem::molarMass(nCompIdx);
Scalar X2 = massFractionW[nCompIdx];
Scalar avgMolarMass = M1*M2/(M2 + X2*(M1 - M2));
// convert mass to mole fractions and set the fluid state
ParentType::fluidState_.setMoleFraction(wPhaseIdx, wCompIdx, massFractionW[wCompIdx]*avgMolarMass/M1);
ParentType::fluidState_.setMoleFraction(wPhaseIdx, nCompIdx, massFractionW[nCompIdx]*avgMolarMass/M2);
// setMassFraction() has only to be called 1-numComponents times
ParentType::fluidState_.setMassFraction(wPhaseIdx, nCompIdx, priVars[switchIdx]);
// TODO give values for non-existing nonwetting phase
Scalar xnH2O = FluidSystem::equilibriumMoleFraction(ParentType::fluidState_, paramCache, nPhaseIdx);
......
......@@ -51,9 +51,13 @@ public:
CompositionalFluidState()
{
// set the composition to 0
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
{
for (int compIdx = 0; compIdx < numComponents; ++compIdx)
{
moleFraction_[phaseIdx][compIdx] = 0;
massFraction_[phaseIdx][compIdx] = 0;
}
averageMolarMass_[phaseIdx] = 0;
sumMoleFractions_[phaseIdx] = 0;
......@@ -87,13 +91,7 @@ public:
* \brief The mass fraction of a component in a phase []
*/
Scalar massFraction(int phaseIdx, int compIdx) const
{
return
std::abs(sumMoleFractions_[phaseIdx])
* moleFraction_[phaseIdx][compIdx]
* FluidSystem::molarMass(compIdx)
/ std::max(1e-40, std::abs(averageMolarMass_[phaseIdx]));
}
{ return massFraction_[phaseIdx][compIdx]; }
/*!
* \brief The average molar mass of a fluid phase [kg/mol]
......@@ -282,6 +280,55 @@ public:
sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
}
for (int compJIdx = 0; compJIdx < numComponents; ++compJIdx)
{
massFraction_[phaseIdx][compJIdx] =
std::abs(sumMoleFractions_[phaseIdx])
* moleFraction_[phaseIdx][compJIdx]
* FluidSystem::molarMass(compJIdx)
/ std::max(1e-40, std::abs(averageMolarMass_[phaseIdx]));
}
}
/*!
* \brief Set the mass fraction of a component in a phase []
* and update the average molar mass [kg/mol] according
* to the current composition of the phase
*/
void setMassFraction(int phaseIdx, int compIdx, Scalar value)
{
Valgrind::CheckDefined(value);
Valgrind::SetDefined(sumMoleFractions_[phaseIdx]);
Valgrind::SetDefined(averageMolarMass_[phaseIdx]);
Valgrind::SetDefined(moleFraction_[phaseIdx][compIdx]);
Valgrind::SetDefined(massFraction_[phaseIdx][compIdx]);
if (numComponents != 2)
DUNE_THROW(Dune::NotImplemented, "This currently only works for 2 components and for setting the mass fractions of the transported component.");
else
{
// calculate average molar mass of the gas phase
Scalar M1 = FluidSystem::molarMass(compIdx);
Scalar M2 = FluidSystem::molarMass(1-compIdx);
Scalar X2 = 1.0-value;
Scalar avgMolarMass = M1*M2/(M2 + X2*(M1 - M2));
massFraction_[phaseIdx][compIdx] = value;
massFraction_[phaseIdx][1-compIdx] = 1.0-value;
moleFraction_[phaseIdx][compIdx] = massFraction_[phaseIdx][compIdx]
* avgMolarMass / M1;
moleFraction_[phaseIdx][1-compIdx] = 1.0-moleFraction_[phaseIdx][compIdx];
// re-calculate the mean molar mass
sumMoleFractions_[phaseIdx] = 0.0;
averageMolarMass_[phaseIdx] = 0.0;
for (int compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
}
}
}
/*!
......@@ -322,6 +369,7 @@ public:
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
for (int compIdx = 0; compIdx < numComponents; ++compIdx) {
Valgrind::CheckDefined(moleFraction_[phaseIdx][compIdx]);
Valgrind::CheckDefined(massFraction_[phaseIdx][compIdx]);
Valgrind::CheckDefined(fugacityCoefficient_[phaseIdx][compIdx]);
}
Valgrind::CheckDefined(averageMolarMass_[phaseIdx]);
......@@ -338,6 +386,7 @@ public:
protected:
Scalar moleFraction_[numPhases][numComponents];
Scalar massFraction_[numPhases][numComponents];
Scalar fugacityCoefficient_[numPhases][numComponents];
Scalar averageMolarMass_[numPhases];
......
......@@ -490,19 +490,8 @@ private:
{
fluidState.setTemperature(refTemperature());
fluidState.setPressure(phaseIdx, refPressure());
Scalar massFraction[numComponents];
massFraction[transportCompIdx] = refMassfrac();
massFraction[phaseCompIdx] = 1 - massFraction[transportCompIdx];
// calculate average molar mass of the gas phase
Scalar M1 = FluidSystem::molarMass(transportCompIdx);
Scalar M2 = FluidSystem::molarMass(phaseCompIdx);
Scalar X2 = massFraction[phaseCompIdx];
Scalar massToMoleDenominator = M2 + X2*(M1 - M2);
fluidState.setMoleFraction(phaseIdx, transportCompIdx, massFraction[transportCompIdx]*M2/massToMoleDenominator);
fluidState.setMoleFraction(phaseIdx, phaseCompIdx, massFraction[phaseCompIdx]*M1/massToMoleDenominator);
// setMassFraction() has only to be called 1-numComponents times
fluidState.setMassFraction(phaseIdx, transportCompIdx, refMassfrac());
}
// can be used for the variation of a boundary condition
......
......@@ -446,9 +446,9 @@ public:
elemVolVarsCur1,
/*onBoundary=*/true);
advectiveWaterVaporFlux += computeAdvectiveVaporFluxes1(elemVolVarsCur1, boundaryVars1, vertInElem1);
diffusiveWaterVaporFlux += computeDiffusiveVaporFluxes1(elemVolVarsCur1, boundaryVars1, vertInElem1);
totalWaterComponentFlux += firstVertexDefect[vertInElem1][transportEqIdx1];
advectiveWaterVaporFlux += computeAdvectiveVaporFluxes1(elemVolVarsCur1, boundaryVars1, vertInElem1);
diffusiveWaterVaporFlux += computeDiffusiveVaporFluxes1(elemVolVarsCur1, boundaryVars1, vertInElem1);
totalWaterComponentFlux += firstVertexDefect[vertInElem1][transportEqIdx1];
}
}
} // end loop over element faces on interface
......
......@@ -462,19 +462,8 @@ private:
{
fluidState.setTemperature(refTemperature());
fluidState.setPressure(phaseIdx, refPressure());
Scalar massFraction[numComponents];
massFraction[transportCompIdx] = refMassfrac();
massFraction[phaseCompIdx] = 1 - massFraction[transportCompIdx];
// calculate average molar mass of the gas phase
Scalar M1 = FluidSystem::molarMass(transportCompIdx);
Scalar M2 = FluidSystem::molarMass(phaseCompIdx);
Scalar X2 = massFraction[phaseCompIdx];
Scalar massToMoleDenominator = M2 + X2*(M1 - M2);
fluidState.setMoleFraction(phaseIdx, transportCompIdx, massFraction[transportCompIdx]*M2/massToMoleDenominator);
fluidState.setMoleFraction(phaseIdx, phaseCompIdx, massFraction[phaseCompIdx]*M1/massToMoleDenominator);
// setMassFraction() has only to be called 1-numComponents times
fluidState.setMassFraction(phaseIdx, transportCompIdx, refMassfrac());
}
// can be used for the variation of a boundary condition
......
......@@ -409,19 +409,8 @@ private:
{
fluidState.setTemperature(refTemperature());
fluidState.setPressure(phaseIdx, refPressure());
Scalar massFraction[numComponents];
massFraction[transportCompIdx] = refMassfrac();
massFraction[phaseCompIdx] = 1 - massFraction[transportCompIdx];
// calculate average molar mass of the gas phase
Scalar M1 = FluidSystem::molarMass(transportCompIdx);
Scalar M2 = FluidSystem::molarMass(phaseCompIdx);
Scalar X2 = massFraction[phaseCompIdx];
Scalar massToMoleDenominator = M2 + X2*(M1 - M2);
fluidState.setMoleFraction(phaseIdx, transportCompIdx, massFraction[transportCompIdx]*M2/massToMoleDenominator);
fluidState.setMoleFraction(phaseIdx, phaseCompIdx, massFraction[phaseCompIdx]*M1/massToMoleDenominator);
// setMassFraction() has only to be called 1-numComponents times
fluidState.setMassFraction(phaseIdx, transportCompIdx, refMassfrac());
}
// can be used for the variation of a boundary condition
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment