diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C index 0d75841d67..7df7108328 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C @@ -407,6 +407,20 @@ Foam::PtrList Foam::phaseSystem::dmdts() const } +bool Foam::phaseSystem::incompressible() const +{ + forAll(phaseModels_, phasei) + { + if (!phaseModels_[phasei].incompressible()) + { + return false; + } + } + + return true; +} + + bool Foam::phaseSystem::implicitPhasePressure(const phaseModel& phase) const { return false; diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H index b048d8d8e5..a90e067ffe 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2015-2019 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2015-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -421,6 +421,9 @@ public: //- Return the mass transfer rates for each phase virtual PtrList dmdts() const; + //- Return incompressibility + bool incompressible() const; + // Transfers diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/createFields.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/createFields.H index 0879cca91b..bca1f2ba5c 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/createFields.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/createFields.H @@ -38,14 +38,30 @@ volScalarField p_rgh label pRefCell = 0; scalar pRefValue = 0.0; -setRefCell -( - p, - p_rgh, - pimple.dict(), - pRefCell, - pRefValue -); +if (fluid.incompressible()) +{ + p = max(p_rgh + fluid.rho()*gh, pMin); + + if (p_rgh.needReference()) + { + setRefCell + ( + p, + p_rgh, + pimple.dict(), + pRefCell, + pRefValue + ); + + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + p_rgh = p - fluid.rho()*gh; + } +} mesh.setFluxRequired(p_rgh.name()); PtrList rAUs; diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H index be9fee7a1c..3fe95ee8b2 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H @@ -335,6 +335,11 @@ while (pimple.correct()) pEqn += pEqnComps[phasei]; } + if (fluid.incompressible()) + { + pEqn.setReference(pRefCell, getRefCellValue(p_rgh, pRefCell)); + } + pEqn.solve(); } @@ -393,6 +398,17 @@ while (pimple.correct()) // Update and limit the static pressure p = max(p_rgh + rho*gh, pMin); + // Account for static pressure reference + if (p_rgh.needReference() && fluid.incompressible()) + { + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + } + // Limit p_rgh p_rgh = p - rho*gh; diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pUf/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pUf/pEqn.H index 8f05c5e1df..4999134331 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pUf/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pUf/pEqn.H @@ -324,6 +324,11 @@ while (pimple.correct()) pEqn += pEqnComps[phasei]; } + if (fluid.incompressible()) + { + pEqn.setReference(pRefCell, getRefCellValue(p_rgh, pRefCell)); + } + pEqn.solve(); } @@ -370,6 +375,17 @@ while (pimple.correct()) // Update and limit the static pressure p = max(p_rgh + rho*gh, pMin); + // Account for static pressure reference + if (p_rgh.needReference() && fluid.incompressible()) + { + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + } + // Limit p_rgh p_rgh = p - rho*gh; diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/createFields.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/createFields.H index a38ec92ee8..69ad38e6ce 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/createFields.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/createFields.H @@ -37,14 +37,30 @@ volScalarField p_rgh label pRefCell = 0; scalar pRefValue = 0.0; -setRefCell -( - p, - p_rgh, - pimple.dict(), - pRefCell, - pRefValue -); +if (fluid.incompressible()) +{ + p = max(p_rgh + fluid.rho()*gh, pMin); + + if (p_rgh.needReference()) + { + setRefCell + ( + p, + p_rgh, + pimple.dict(), + pRefCell, + pRefValue + ); + + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + p_rgh = p - fluid.rho()*gh; + } +} mesh.setFluxRequired(p_rgh.name()); PtrList rAUs; diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/pEqn.H index 4e4ed097e3..0d00e8a16e 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/pEqn.H @@ -286,6 +286,11 @@ while (pimple.correct()) - fvm::laplacian(rAUf, p_rgh) ); + if (fluid.incompressible()) + { + pEqnIncomp.setReference(pRefCell, getRefCellValue(p_rgh, pRefCell)); + } + // Solve solve(pEqnIncomp + pEqnComp1 + pEqnComp2); @@ -365,6 +370,17 @@ while (pimple.correct()) // Update and limit the static pressure p = max(p_rgh + rho*gh, pMin); + // Account for static pressure reference + if (p_rgh.needReference() && fluid.incompressible()) + { + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + } + // Limit p_rgh p_rgh = p - rho*gh; diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/pEqn.H index f88e7af343..e679cbd107 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/pEqn.H @@ -298,6 +298,11 @@ while (pimple.correct()) - fvm::laplacian(rAUf, p_rgh) ); + if (fluid.incompressible()) + { + pEqnIncomp.setReference(pRefCell, getRefCellValue(p_rgh, pRefCell)); + } + // Solve solve(pEqnIncomp + pEqnComp1 + pEqnComp2); @@ -347,6 +352,17 @@ while (pimple.correct()) // Update and limit the static pressure p = max(p_rgh + rho*gh, pMin); + // Account for static pressure reference + if (p_rgh.needReference() && fluid.incompressible()) + { + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + } + // Limit p_rgh p_rgh = p - rho*gh; diff --git a/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/trickleBed/system/fvSolution b/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/trickleBed/system/fvSolution index a3cac2f345..bb274db11e 100644 --- a/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/trickleBed/system/fvSolution +++ b/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/trickleBed/system/fvSolution @@ -61,8 +61,6 @@ PIMPLE nOuterCorrectors 3; nCorrectors 1; nNonOrthogonalCorrectors 0; - pRefCell 0; - pRefValue 1e5; partialElimination true; }