diff --git a/dumux/porousmediumflow/2pnc/implicit/primaryvariableswitch.hh b/dumux/porousmediumflow/2pnc/implicit/primaryvariableswitch.hh index 145fb41b99777f76b0529895c57475f41f7f8d4c..941767ba9601d3bccc8a4767ec87c44b3dba1d66 100644 --- a/dumux/porousmediumflow/2pnc/implicit/primaryvariableswitch.hh +++ b/dumux/porousmediumflow/2pnc/implicit/primaryvariableswitch.hh @@ -95,92 +95,92 @@ protected: if (this->wasSwitched_[dofIdxGlobal]) Smin = -0.01; - //if saturation of liquid phase is smaller 0 switch + //if saturation of wetting phase is smaller 0 switch if (volVars.saturation(wPhaseIdx) <= Smin) { wouldSwitch = true; - //liquid phase has to disappear - std::cout << "Liquid Phase disappears at vertex " << dofIdxGlobal - << ", coordinated: " << globalPos << ", Sl: " + //wetting phase has to disappear + std::cout << "Wetting Phase disappears at vertex " << dofIdxGlobal + << ", coordinated: " << globalPos << ", Sw: " << volVars.saturation(wPhaseIdx) << std::endl; newPhasePresence = nPhaseOnly; //switch not depending on formulation - //switch "Sl" to "xgH20" + //switch "Sw" to "xn20" priVars[switchIdx] = volVars.moleFraction(nPhaseIdx, wCompIdx /*H2O*/); - //switch all secondary components to mole fraction in gas phase + //switch all secondary components to mole fraction in non-wetting phase for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx) priVars[compIdx] = volVars.moleFraction(nPhaseIdx,compIdx); } - //if saturation of gas phase is smaller than 0 switch + //if saturation of non-wetting phase is smaller than 0 switch else if (volVars.saturation(nPhaseIdx) <= Smin) { wouldSwitch = true; - //gas phase has to disappear - std::cout << "Gas Phase disappears at vertex " << dofIdxGlobal - << ", coordinated: " << globalPos << ", Sg: " + //non-wetting phase has to disappear + std::cout << "Non-wetting Phase disappears at vertex " << dofIdxGlobal + << ", coordinated: " << globalPos << ", Sn: " << volVars.saturation(nPhaseIdx) << std::endl; newPhasePresence = wPhaseOnly; - //switch "Sl" to "xlN2" + //switch "Sn" to "xwN2" priVars[switchIdx] = volVars.moleFraction(wPhaseIdx, nCompIdx /*N2*/); } } else if (phasePresence == nPhaseOnly) { - Scalar xlmax = 1; - Scalar sumxl = 0; - //Calculate sum of mole fractions in the hypothetical liquid phase + Scalar xwmax = 1; + Scalar sumxw = 0; + //Calculate sum of mole fractions in the hypothetical wetting phase for (int compIdx = 0; compIdx < numComponents; compIdx++) { - sumxl += volVars.moleFraction(wPhaseIdx, compIdx); + sumxw += volVars.moleFraction(wPhaseIdx, compIdx); } - if (sumxl > xlmax) + if (sumxw > xwmax) wouldSwitch = true; if (this->wasSwitched_[dofIdxGlobal]) - xlmax *=1.02; - //liquid phase appears if sum is larger than one - if (sumxl/*sum of mole fractions*/ > xlmax/*1*/) + xwmax *=1.02; + //wetting phase appears if sum is larger than one + if (sumxw/*sum of mole fractions*/ > xwmax/*1*/) { - std::cout << "Liquid Phase appears at vertex " << dofIdxGlobal - << ", coordinated: " << globalPos << ", sumxl: " - << sumxl << std::endl; + std::cout << "Wetting Phase appears at vertex " << dofIdxGlobal + << ", coordinated: " << globalPos << ", sumxw: " + << sumxw << std::endl; newPhasePresence = bothPhases; - //saturation of the liquid phase set to 0.0001 (if formulation pgSl and vice versa) + //saturation of the wetting phase set to 0.0001 (if formulation pnsw and vice versa) if (formulation == pnsw) priVars[switchIdx] = 0.0001; else if (formulation == pwsn) priVars[switchIdx] = 0.9999; - //switch all secondary components back to liquid mole fraction + //switch all secondary components back to wetting mole fraction for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx) priVars[compIdx] = volVars.moleFraction(wPhaseIdx,compIdx); } } else if (phasePresence == wPhaseOnly) { - Scalar xgmax = 1; - Scalar sumxg = 0; - //Calculate sum of mole fractions in the hypothetical liquid phase + Scalar xnmax = 1; + Scalar sumxn = 0; + //Calculate sum of mole fractions in the hypothetical wetting phase for (int compIdx = 0; compIdx < numComponents; compIdx++) { - sumxg += volVars.moleFraction(nPhaseIdx, compIdx); + sumxn += volVars.moleFraction(nPhaseIdx, compIdx); } - if (sumxg > xgmax) + if (sumxn > xnmax) wouldSwitch = true; if (this->wasSwitched_[dofIdxGlobal]) - xgmax *=1.02; - //liquid phase appears if sum is larger than one - if (sumxg > xgmax) + xnmax *=1.02; + //wetting phase appears if sum is larger than one + if (sumxn > xnmax) { - std::cout << "Gas Phase appears at vertex " << dofIdxGlobal - << ", coordinated: " << globalPos << ", sumxg: " - << sumxg << std::endl; + std::cout << "Non-wetting Phase appears at vertex " << dofIdxGlobal + << ", coordinated: " << globalPos << ", sumxn: " + << sumxn << std::endl; newPhasePresence = bothPhases; - //saturation of the liquid phase set to 0.9999 (if formulation pgSl and vice versa) + //saturation of the wetting phase set to 0.9999 (if formulation pnsw and vice versa) if (formulation == pnsw) priVars[switchIdx] = 0.9999; else if (formulation == pwsn) diff --git a/dumux/porousmediumflow/2pnc/implicit/volumevariables.hh b/dumux/porousmediumflow/2pnc/implicit/volumevariables.hh index 8783a1524001eee81d69ed1c790ecdc5d0b62dfc..01b56f98cb79c729f73efe2c3ce5359c3024ca3b 100644 --- a/dumux/porousmediumflow/2pnc/implicit/volumevariables.hh +++ b/dumux/porousmediumflow/2pnc/implicit/volumevariables.hh @@ -186,21 +186,21 @@ public: // set the saturations ///////////// - Scalar Sg; + Scalar Sn; if (phasePresence == nPhaseOnly) { - Sg = 1.0; + Sn = 1.0; } else if (phasePresence == wPhaseOnly) { - Sg = 0.0; + Sn = 0.0; } else if (phasePresence == bothPhases) { if (formulation == pwsn) - Sg = priVars[switchIdx]; + Sn = priVars[switchIdx]; else if (formulation == pnsw) - Sg = 1.0 - priVars[switchIdx]; + Sn = 1.0 - priVars[switchIdx]; else DUNE_THROW(Dune::InvalidStateException, "Formulation: " << formulation << " is invalid."); } @@ -209,8 +209,8 @@ public: DUNE_THROW(Dune::InvalidStateException, "phasePresence: " << phasePresence << " is invalid."); } - fluidState.setSaturation(nPhaseIdx, Sg); - fluidState.setSaturation(wPhaseIdx, 1.0 - Sg); + fluidState.setSaturation(nPhaseIdx, Sn); + fluidState.setSaturation(wPhaseIdx, 1.0 - Sn); ///////////// // set the pressures of the fluid phases @@ -218,7 +218,7 @@ public: // calculate capillary pressure const auto& materialParams = problem.spatialParams().materialLawParams(element, scv, elemSol); - auto pc = MaterialLaw::pc(materialParams, 1 - Sg); + auto pc = MaterialLaw::pc(materialParams, 1 - Sn); // extract the pressures if (formulation == pwsn) @@ -234,7 +234,7 @@ public: // Here we check for (p_g - pc) in order to ensure that (p_l > 0) if (priVars[pressureIdx] - pc < 0.0) { - std::cout<< "p_g: "<< priVars[pressureIdx]<<" Cap_press: "<< pc << std::endl; + std::cout<< "p_n: "<< priVars[pressureIdx]<<" Cap_press: "<< pc << std::endl; DUNE_THROW(NumericalProblem,"Capillary pressure is too high"); } fluidState.setPressure(wPhaseIdx, priVars[pressureIdx] - pc); @@ -313,19 +313,20 @@ public: { // only the wetting phase is present, i.e. wetting phase // composition is stored explicitly. - // extract _mass_ fractions in the nonwetting phase + // extract _mass_ fractions in the non-wetting phase Dune::FieldVector<Scalar, numComponents> moleFrac; + moleFrac[nCompIdx] = priVars[switchIdx]; + for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx) moleFrac[compIdx] = priVars[compIdx]; - moleFrac[nCompIdx] = priVars[switchIdx]; - Scalar sumMoleFracNotWater = 0; + Scalar sumMoleFracOtherComponents = 0; for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx) - sumMoleFracNotWater+=moleFrac[compIdx]; + sumMoleFracOtherComponents += moleFrac[compIdx]; - sumMoleFracNotWater += moleFrac[nCompIdx]; - moleFrac[wCompIdx] = 1 -sumMoleFracNotWater; + sumMoleFracOtherComponents += moleFrac[nCompIdx]; + moleFrac[wCompIdx] = 1 - sumMoleFracOtherComponents; // convert mass to mole fractions and set the fluid state for (int compIdx=0; compIdx<numComponents; ++compIdx) diff --git a/test/porousmediumflow/2pnc/implicit/fuelcellproblem.hh b/test/porousmediumflow/2pnc/implicit/fuelcellproblem.hh index e854362a31f926ba771ce8c5f98294d8b043ddce..05e560b2b3c74239d8f99f3f921125e65e40d536 100644 --- a/test/porousmediumflow/2pnc/implicit/fuelcellproblem.hh +++ b/test/porousmediumflow/2pnc/implicit/fuelcellproblem.hh @@ -245,9 +245,9 @@ public: if(onUpperBoundary_(globalPos)) { - Scalar pg = 1.0e5; - priVars[pressureIdx] = pg; - priVars[switchIdx] = 0.3;//Sl for bothPhases + Scalar pn = 1.0e5; + priVars[pressureIdx] = pn; + priVars[switchIdx] = 0.3;//Sw for bothPhases priVars[switchIdx+1] = pO2Inlet_/4.315e9; //moleFraction xlO2 for bothPhases } @@ -346,9 +346,9 @@ private: PrimaryVariables priVars(0.0); priVars.setState(Indices::bothPhases); - Scalar pg = 1.0e5; - priVars[pressureIdx] = pg; - priVars[switchIdx] = 0.3;//Sl for bothPhases + Scalar pn = 1.0e5; + priVars[pressureIdx] = pn; + priVars[switchIdx] = 0.3;//Sw for bothPhases priVars[switchIdx+1] = pO2Inlet_/4.315e9; //moleFraction xlO2 for bothPhases return priVars;