From 0196822968ef777c990cffa86d6b5b54a6ff16ba Mon Sep 17 00:00:00 2001 From: Henry Weller <http://cfd.direct> Date: Mon, 15 Feb 2016 18:00:20 +0000 Subject: [PATCH] MPPICFoam PackingModels: Implicit: Added limiting to the calculation of the correction flux. Vastly reduces the scattering and churning behaviour of packed beds. Development provided by Will Bainbridge <github.com/will-bainbridge> See also http://www.openfoam.org/mantisbt/view.php?id=1994 --- .../MPPIC/PackingModels/Implicit/Implicit.C | 102 +++++++++++++++--- .../MPPIC/PackingModels/Implicit/Implicit.H | 5 +- .../constant/kinematicCloudProperties | 1 + .../column/constant/kinematicCloudProperties | 1 + .../cyclone/constant/kinematicCloudProperties | 1 + 5 files changed, 93 insertions(+), 17 deletions(-) diff --git a/src/lagrangian/intermediate/submodels/MPPIC/PackingModels/Implicit/Implicit.C b/src/lagrangian/intermediate/submodels/MPPIC/PackingModels/Implicit/Implicit.C index 5cb2092818c..ab2b4579915 100644 --- a/src/lagrangian/intermediate/submodels/MPPIC/PackingModels/Implicit/Implicit.C +++ b/src/lagrangian/intermediate/submodels/MPPIC/PackingModels/Implicit/Implicit.C @@ -48,6 +48,7 @@ Foam::PackingModels::Implicit<CloudType>::Implicit ), phiCorrect_(NULL), uCorrect_(NULL), + applyLimiting_(this->coeffDict().lookup("applyLimiting")), applyGravity_(this->coeffDict().lookup("applyGravity")), alphaMin_(readScalar(this->coeffDict().lookup("alphaMin"))), rhoMin_(readScalar(this->coeffDict().lookup("rhoMin"))) @@ -66,6 +67,7 @@ Foam::PackingModels::Implicit<CloudType>::Implicit alpha_(cm.alpha_), phiCorrect_(cm.phiCorrect_()), uCorrect_(cm.uCorrect_()), + applyLimiting_(cm.applyLimiting_), applyGravity_(cm.applyGravity_), alphaMin_(cm.alphaMin_), rhoMin_(cm.rhoMin_) @@ -102,6 +104,11 @@ void Foam::PackingModels::Implicit<CloudType>::cacheFields(const bool store) ( cloudName + ":rhoAverage" ); + const AveragingMethod<vector>& uAverage = + mesh.lookupObject<AveragingMethod<vector> > + ( + cloudName + ":uAverage" + ); const AveragingMethod<scalar>& uSqrAverage = mesh.lookupObject<AveragingMethod<scalar>> ( @@ -135,13 +142,6 @@ void Foam::PackingModels::Implicit<CloudType>::cacheFields(const bool store) rho.internalField() = max(rhoAverage.internalField(), rhoMin_); rho.correctBoundaryConditions(); - //Info << " x: " << mesh.C().internalField().component(2) << endl; - //Info << " alpha: " << alpha_.internalField() << endl; - //Info << "alphaOld: " << alpha_.oldTime().internalField() << endl; - //Info << " rho: " << rho.internalField() << endl; - //Info << endl; - - // Stress field // ~~~~~~~~~~~~ @@ -172,6 +172,24 @@ void Foam::PackingModels::Implicit<CloudType>::cacheFields(const bool store) tauPrime.correctBoundaryConditions(); + // Gravity flux + // ~~~~~~~~~~~~ + + tmp<surfaceScalarField> phiGByA; + + if (applyGravity_) + ( + phiGByA = tmp<surfaceScalarField> + ( + new surfaceScalarField + ( + "phiGByA", + deltaT*(g & mesh.Sf())*fvc::interpolate(1.0 - rhoc/rho) + ) + ) + ); + + // Implicit solution for the volume fraction // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -191,14 +209,7 @@ void Foam::PackingModels::Implicit<CloudType>::cacheFields(const bool store) if (applyGravity_) { - surfaceScalarField - phiGByA - ( - "phiGByA", - deltaT*(g & mesh.Sf())*fvc::interpolate(1.0 - rhoc/rho) - ); - - alphaEqn += fvm::div(phiGByA, alpha_); + alphaEqn += fvm::div(phiGByA(), alpha_); } alphaEqn.solve(); @@ -217,6 +228,67 @@ void Foam::PackingModels::Implicit<CloudType>::cacheFields(const bool store) ) ); + // limit the correction flux + if (applyLimiting_) + { + volVectorField U + ( + IOobject + ( + cloudName + ":U", + this->owner().db().time().timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh, + dimensionedVector("zero", dimVelocity, vector::zero), + fixedValueFvPatchField<vector>::typeName + ); + U.internalField() = uAverage.internalField(); + U.correctBoundaryConditions(); + + surfaceScalarField phi + ( + cloudName + ":phi", + linearInterpolate(U) & mesh.Sf() + ); + + if (applyGravity_) + { + phiCorrect_() -= phiGByA(); + } + + forAll(phiCorrect_(), faceI) + { + // Current and correction fluxes + const scalar phiCurr = phi[faceI]; + scalar& phiCorr = phiCorrect_()[faceI]; + + // Don't limit if the correction is in the opposite direction to + // the flux. We need all the help we can get in this state. + if (phiCurr*phiCorr < 0) + {} + + // If the correction and the flux are in the same direction then + // don't apply any more correction than is already present in + // the flux. + else if (phiCorr > 0) + { + phiCorr = max(phiCorr - phiCurr, 0); + } + else + { + phiCorr = min(phiCorr - phiCurr, 0); + } + } + + if (applyGravity_) + { + phiCorrect_() += phiGByA(); + } + } + // correction velocity uCorrect_ = tmp<volVectorField> ( diff --git a/src/lagrangian/intermediate/submodels/MPPIC/PackingModels/Implicit/Implicit.H b/src/lagrangian/intermediate/submodels/MPPIC/PackingModels/Implicit/Implicit.H index 9f90dabcbf5..dc311700040 100644 --- a/src/lagrangian/intermediate/submodels/MPPIC/PackingModels/Implicit/Implicit.H +++ b/src/lagrangian/intermediate/submodels/MPPIC/PackingModels/Implicit/Implicit.H @@ -58,8 +58,6 @@ class Implicit : public PackingModel<CloudType> { -private: - //- Private data //- Volume fraction field @@ -71,6 +69,9 @@ private: //- Correction cell-centred velocity tmp<volVectorField> uCorrect_; + //- Flag to indicate whether implicit limiting is applied + Switch applyLimiting_; + //- Flag to indicate whether gravity is applied Switch applyGravity_; diff --git a/tutorials/lagrangian/MPPICFoam/Goldschmidt/constant/kinematicCloudProperties b/tutorials/lagrangian/MPPICFoam/Goldschmidt/constant/kinematicCloudProperties index fdb104bc2d3..9f46325f606 100644 --- a/tutorials/lagrangian/MPPICFoam/Goldschmidt/constant/kinematicCloudProperties +++ b/tutorials/lagrangian/MPPICFoam/Goldschmidt/constant/kinematicCloudProperties @@ -148,6 +148,7 @@ subModels { alphaMin 0.0001; rhoMin 1.0; + applyLimiting true; applyGravity false; particleStressModel { diff --git a/tutorials/lagrangian/MPPICFoam/column/constant/kinematicCloudProperties b/tutorials/lagrangian/MPPICFoam/column/constant/kinematicCloudProperties index f086e389ec1..00d205137b5 100644 --- a/tutorials/lagrangian/MPPICFoam/column/constant/kinematicCloudProperties +++ b/tutorials/lagrangian/MPPICFoam/column/constant/kinematicCloudProperties @@ -140,6 +140,7 @@ subModels { alphaMin 0.0001; rhoMin 1.0; + applyLimiting true; applyGravity true; particleStressModel { diff --git a/tutorials/lagrangian/MPPICFoam/cyclone/constant/kinematicCloudProperties b/tutorials/lagrangian/MPPICFoam/cyclone/constant/kinematicCloudProperties index 78bac15385f..d72fd631630 100644 --- a/tutorials/lagrangian/MPPICFoam/cyclone/constant/kinematicCloudProperties +++ b/tutorials/lagrangian/MPPICFoam/cyclone/constant/kinematicCloudProperties @@ -145,6 +145,7 @@ subModels { alphaMin 0.0001; rhoMin 1.0; + applyLimiting true; applyGravity false; particleStressModel { -- GitLab