From a7f747a4abe6be61cf519863fa163a3e413b8f33 Mon Sep 17 00:00:00 2001 From: tlichtenegger Date: Mon, 9 Sep 2019 11:42:45 +0200 Subject: [PATCH] Added stickiness to dust model; also some cleaning up. --- .../subModels/forceModel/Fines/FinesFields.C | 153 ++++++++++++------ .../subModels/forceModel/Fines/FinesFields.H | 16 +- 2 files changed, 116 insertions(+), 53 deletions(-) diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.C b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.C index 31b5e008..7bb74c3f 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.C +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.C @@ -47,11 +47,13 @@ FinesFields::FinesFields verbose_(propsDict_.lookupOrDefault("verbose",false)), clogKin_(propsDict_.lookupOrDefault("kineticClogging",false)), clogStick_(propsDict_.lookupOrDefault("stickyClogging",false)), + movingBed_(propsDict_.lookupOrDefault("movingBed",true)), + fluxFieldName_(propsDict_.lookupOrDefault("fluxFieldName","phi")), + phi_(sm.mesh().lookupObject (fluxFieldName_)), velFieldName_(propsDict_.lookupOrDefault("velFieldName","U")), U_(sm.mesh().lookupObject (velFieldName_)), voidfractionFieldName_(propsDict_.lookupOrDefault("voidfractionFieldName","voidfraction")), - // this is probably really bad - voidfraction_(const_cast(sm.mesh().lookupObject (voidfractionFieldName_))), + voidfraction_(sm.mesh().lookupObject (voidfractionFieldName_)), UsFieldName_(propsDict_.lookupOrDefault("granVelFieldName","Us")), UsField_(sm.mesh().lookupObject (UsFieldName_)), pFieldName_(propsDict_.lookupOrDefault("pFieldName","p")), @@ -202,7 +204,6 @@ FinesFields::FinesFields ), sm.mesh(), dimensionedScalar("zero", dimensionSet(1,0,-1,0,0), 0) - //dimensionedVector("zero", dimensionSet(1,-2,-1,0,0), vector::zero) ), uDyn_ ( IOobject @@ -221,8 +222,10 @@ FinesFields::FinesFields rhoFine_("rhoFine",dimensionSet(1,-3,0,0,0),0.0), g_("g",dimensionSet(0,1,-2,0,0),vector(0,0,-9.81)), alphaDynMax_(0.1), - alphaMax_(readScalar(propsDict_.lookup ("alphaMax"))), - critVoidfraction_(readScalar(propsDict_.lookup ("critVoidfraction"))), + alphaMax_(propsDict_.lookupOrDefault("alphaMax",0.95)), + alphaMinClog_(propsDict_.lookupOrDefault("alphaMinClog",0.3)), + critVoidfraction_(propsDict_.lookupOrDefault("critVoidfraction", 0.05)), + deltaT_(voidfraction_.mesh().time().deltaTValue()), depositionLength_(readScalar(propsDict_.lookup ("depositionLength"))), exponent_(-1.33), nCrit_(0.0), @@ -281,26 +284,48 @@ FinesFields::FinesFields if (clogStick_) { uBind_ = readScalar(propsDict_.lookup ("uBind")); + + uReconstructed_.set + ( + new volVectorField + ( + IOobject + ( + "uReconstructed", + sm.mesh().time().timeName(), + sm.mesh(), + IOobject::NO_READ, + IOobject::NO_WRITE + ), + U_//sm.mesh() + ) + ); + } - if(verbose_) + if (verbose_) { - alphaG_.writeOpt() = IOobject::AUTO_WRITE; - alphaG_.write(); - alphaP_.writeOpt() = IOobject::AUTO_WRITE; - alphaP_.write(); - deltaAlpha_.writeOpt() = IOobject::AUTO_WRITE; - deltaAlpha_.write(); - dHydMix_.writeOpt() = IOobject::AUTO_WRITE; - dHydMix_.write(); - DragCoeff_.writeOpt() = IOobject::AUTO_WRITE; - DragCoeff_.write(); - dSauterMix_.writeOpt() = IOobject::AUTO_WRITE; - dSauterMix_.write(); - FanningCoeff_.writeOpt() = IOobject::AUTO_WRITE; - FanningCoeff_.write(); - Froude_.writeOpt() = IOobject::AUTO_WRITE; - Froude_.write(); + alphaG_.writeOpt() = IOobject::AUTO_WRITE; + alphaG_.write(); + alphaP_.writeOpt() = IOobject::AUTO_WRITE; + alphaP_.write(); + deltaAlpha_.writeOpt() = IOobject::AUTO_WRITE; + deltaAlpha_.write(); + dHydMix_.writeOpt() = IOobject::AUTO_WRITE; + dHydMix_.write(); + DragCoeff_.writeOpt() = IOobject::AUTO_WRITE; + DragCoeff_.write(); + dSauterMix_.writeOpt() = IOobject::AUTO_WRITE; + dSauterMix_.write(); + FanningCoeff_.writeOpt() = IOobject::AUTO_WRITE; + FanningCoeff_.write(); + Froude_.writeOpt() = IOobject::AUTO_WRITE; + Froude_.write(); + if (clogStick_) + { + uReconstructed_().writeOpt() = IOobject::AUTO_WRITE; + uReconstructed_().write(); + } } } @@ -315,33 +340,38 @@ FinesFields::~FinesFields() void FinesFields::update() { - if(verbose_) Info << "FinesFields: Updating alphaP.\n" << endl; + if (verbose_) Info << "FinesFields: Updating alphaP.\n" << endl; updateAlphaP(); - if(verbose_) Info << "FinesFields: Updating alphaG.\n" << endl; + if (verbose_) Info << "FinesFields: Updating alphaG.\n" << endl; updateAlphaG(); - if(verbose_) Info << "FinesFields: Calculating source terms.\n" << endl; + if (clogStick_) + { + if (verbose_) Info << "FinesFields: Updating uReconstructed.\n" << endl; + updateUReconstructed(); + } + if (verbose_) Info << "FinesFields: Calculating source terms.\n" << endl; calcSource(); - if(verbose_) Info << "FinesFields: Updating dSauter.\n" << endl; + if (verbose_) Info << "FinesFields: Updating dSauter.\n" << endl; updateDSauter(); - if(verbose_) Info << "FinesFields: Updating dHydMix.\n" << endl; + if (verbose_) Info << "FinesFields: Updating dHydMix.\n" << endl; updateDHydMix(); - if(verbose_) Info << "FinesFields: Updating Froude.\n" << endl; + if (verbose_) Info << "FinesFields: Updating Froude.\n" << endl; updateFroude(); - if(verbose_) Info << "FinesFields: Updating FanningCoeff.\n" << endl; + if (verbose_) Info << "FinesFields: Updating FanningCoeff.\n" << endl; updateFanningCoeff(); - if(verbose_) Info << "FinesFields: Updating DragCoeff.\n" << endl; + if (verbose_) Info << "FinesFields: Updating DragCoeff.\n" << endl; updateDragCoeff(); - if(verbose_) Info << "FinesFields: Updating uDyn.\n" << endl; + if (verbose_) Info << "FinesFields: Updating uDyn.\n" << endl; updateUDyn(); - if(verbose_) Info << "FinesFields: Integrating alphas.\n" << endl; + if (verbose_) Info << "FinesFields: Integrating alphas.\n" << endl; integrateFields(); - if(verbose_) Info << "FinesFields: Update finished.\n" << endl; + if (verbose_) Info << "FinesFields: Update finished.\n" << endl; } void FinesFields::calcSource() { - if(!clogKin_ & !clogStick_) return; + if (!clogKin_ & !clogStick_) return; Sds_.primitiveFieldRef() = 0; deltaAlpha_.primitiveFieldRef() = 0.0; scalar fKin = 0.0; @@ -350,11 +380,14 @@ void FinesFields::calcSource() scalar dmean = 0.0; scalar d1 = 0.0; scalar d2 = 0.0; + scalar magU = 0.0; scalar tauDeposition = 0.0; forAll(Sds_,cellI) { - if(clogKin_) + fKin = 0.0; + fStick = 0.0; + if (clogKin_ && alphaP_[cellI] > alphaMinClog_) // no kinetic cloggig in dilute regions { // calculate everything in units auf dSauter critpore = nCrit_*dFine_.value()/dSauter_[cellI]; @@ -370,9 +403,11 @@ void FinesFields::calcSource() else if (fKin > 1.0) fKin = 1.0; } - if(clogStick_) + if (clogStick_) { - fStick = 1.0 / ( 1.0 + mag(U_[cellI])/uBind_); + magU = mag(uReconstructed_()[cellI]); // use U reconstructed from phi to suppress oscillations at interfaces + fStick = 1.0 / ( 1.0 + magU/uBind_) * alphaP_[cellI] / 0.65; + if (fStick > 1.0) fStick = 1.0; } // at this point, voidfraction is still calculated from the true particle sizes @@ -386,7 +421,8 @@ void FinesFields::calcSource() else { tauDeposition = depositionLength_ / max(mag(U_[cellI]),uMin_); - Sds_[cellI] = deltaAlpha_[cellI] / tauDeposition; + if (tauDeposition < 4.0*deltaT_) tauDeposition = 4.0*deltaT_; // do not deposit all dyn hold-up within less than 3 time steps + Sds_[cellI] = alphaDyn_[cellI] / tauDeposition; } } } @@ -400,11 +436,12 @@ void FinesFields::integrateFields() fvScalarMatrix alphaStEqn ( - fvm::ddt(alphaSt_) - + fvm::div(phiSt,alphaSt_) + fvm::ddt(alphaSt_) == Sds_ ); + if (movingBed_) alphaStEqn += fvm::div(phiSt,alphaSt_); + fvScalarMatrix alphaDynEqn ( fvm::ddt(alphaDyn_) @@ -413,10 +450,11 @@ void FinesFields::integrateFields() == -Sds_ ); + alphaStEqn.solve(); alphaDynEqn.solve(); - if(smoothing_) + if (smoothing_) particleCloud_.smoothingM().smoothen(alphaDyn_); // limit hold-ups, should be done more elegantly @@ -475,7 +513,7 @@ void FinesFields::updateDHydMix() forAll(dHydMix_,cellI) { scalar aPSt = alphaP_[cellI] + alphaSt_[cellI]; - if(aPSt < SMALL || aPSt > 1 - SMALL) + if (aPSt < SMALL || aPSt > 1 - SMALL) dHydMix_[cellI] = SMALL; else dHydMix_[cellI] = 2*(1 - aPSt) / (3*aPSt ) * dSauterMix_[cellI]; @@ -495,11 +533,11 @@ void FinesFields::updateDragCoeff() forAll(DragCoeff_,cellI) { Ref1 = Ref[cellI]; - if(Ref1 <= SMALL) + if (Ref1 <= SMALL) Cd = 24.0 / SMALL; - else if(Ref1 <= 1.0) + else if (Ref1 <= 1.0) Cd = 24.0 / Ref1; - else if(Ref1 <= 1000) + else if (Ref1 <= 1000) Cd = 24 * (1.0 + 0.15 * Foam::pow(Ref1,0.687) ) / Ref1; else Cd = 0.44; @@ -511,11 +549,11 @@ void FinesFields::updateDragCoeff() forAll(DragCoeff_.boundaryField()[patchI], faceI) { Ref1 = Ref.boundaryField()[patchI][faceI]; - if(Ref1 <= SMALL) + if (Ref1 <= SMALL) Cd = 24.0 / SMALL; - else if(Ref1 <= 1.0) + else if (Ref1 <= 1.0) Cd = 24.0 / Ref1; - else if(Ref1 <= 1000) + else if (Ref1 <= 1000) Cd = 24 * (1.0 + 0.15 * Foam::pow(Ref1,0.687) ) / Ref1; else Cd = 0.44; @@ -532,9 +570,9 @@ void FinesFields::updateDSauter() { scalar aP = alphaP_[cellI]; scalar aSt = alphaSt_[cellI]; - if(aSt < SMALL) + if (aSt < SMALL) dSauterMix_[cellI] = dSauter_[cellI]; - else if(aP < SMALL) + else if (aP < SMALL) dSauterMix_[cellI] = dFine_.value(); else dSauterMix_[cellI] = (aP + aSt) / (aP / dSauter_[cellI] + aSt / dFine_.value() ); @@ -577,7 +615,7 @@ void FinesFields::updateUDyn() { scalar mU(mag(U_[cellI])); scalar muDyn(mag(uDyn_[cellI])); - if(muDyn > mU && muDyn > SMALL) + if (muDyn > mU && muDyn > SMALL) { uDyn_[cellI] *= mU / muDyn; } @@ -588,13 +626,24 @@ void FinesFields::updateUDyn() { scalar mU(mag(U_.boundaryField()[patchI][faceI])); scalar muDyn(mag(uDyn_.boundaryField()[patchI][faceI])); - if(muDyn > mU && muDyn > SMALL) + if (muDyn > mU && muDyn > SMALL) { uDyn_.boundaryFieldRef()[patchI][faceI] *= mU / muDyn; } } } +void FinesFields::updateUReconstructed() +{ + if (phi_.dimensions() == dimensionSet(1, 0, -1, 0, 0)) // compressible + { + uReconstructed_() = fvc::reconstruct(phi_) / (rhoG_ * voidfraction_); + } + else + { + uReconstructed_() = fvc::reconstruct(phi_) / voidfraction_; + } +} // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace Foam diff --git a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.H b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.H index d67e79f0..ecc26a77 100644 --- a/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.H +++ b/src/lagrangian/cfdemParticle/subModels/forceModel/Fines/FinesFields.H @@ -50,13 +50,19 @@ private: bool clogStick_; + bool movingBed_; + + word fluxFieldName_; + + const surfaceScalarField& phi_; + word velFieldName_; const volVectorField& U_; word voidfractionFieldName_; - volScalarField& voidfraction_; + const volScalarField& voidfraction_; word UsFieldName_; @@ -98,6 +104,8 @@ private: volVectorField uDyn_; + autoPtr uReconstructed_; // velocity field can show oscillations at porosity steps + dimensionedScalar dFine_; dimensionedScalar diffCoeff_; @@ -112,8 +120,12 @@ private: scalar alphaMax_; + scalar alphaMinClog_; + scalar critVoidfraction_; + scalar deltaT_; + scalar depositionLength_; scalar exponent_; @@ -152,6 +164,8 @@ private: void updateUDyn(); + void updateUReconstructed(); + public: //- Runtime type information