From 29cea780e14eab2c49fe6d53ba2e358284b15c1f Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Fri, 18 Sep 2015 18:55:21 +0100 Subject: [PATCH] reactingMultiphaseEulerFoam: Optimize the handling of optional forces --- .../MomentumTransferPhaseSystem.C | 75 ++++++++++++------- .../MomentumTransferPhaseSystem.H | 14 ++++ .../reactingMultiphaseEulerFoam/pU/pEqn.H | 64 ++++++++++++---- .../bubbleColumn/constant/phaseProperties | 9 +-- 4 files changed, 112 insertions(+), 50 deletions(-) diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C index f43a3a0354..f50bd61f9e 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C @@ -511,17 +511,12 @@ Foam::MomentumTransferPhaseSystem::momentumTransfer() const template -Foam::autoPtr > -Foam::MomentumTransferPhaseSystem::Fs() const +Foam::volVectorField& Foam::MomentumTransferPhaseSystem::setF +( + PtrList& Fs, const label phasei +) const { - autoPtr > tFs - ( - new PtrList(this->phases().size()) - ); - - PtrList& Fs = tFs(); - - forAll(Fs, phasei) + if (!Fs.set(phasei)) { Fs.set ( @@ -543,6 +538,20 @@ Foam::MomentumTransferPhaseSystem::Fs() const ); } + return Fs[phasei]; +} + + +template +Foam::autoPtr > +Foam::MomentumTransferPhaseSystem::Fs() const +{ + autoPtr > tFs + ( + new PtrList(this->phases().size()) + ); + PtrList& Fs = tFs(); + // Add the lift force forAllConstIter ( @@ -555,8 +564,8 @@ Foam::MomentumTransferPhaseSystem::Fs() const const phasePair& pair(this->phasePairs_[liftModelIter.key()]); - Fs[pair.phase1().index()] += F; - Fs[pair.phase2().index()] -= F; + setF(Fs, pair.phase1().index()) += F; + setF(Fs, pair.phase2().index()) -= F; } // Add the wall lubrication force @@ -572,8 +581,8 @@ Foam::MomentumTransferPhaseSystem::Fs() const const phasePair& pair(this->phasePairs_[wallLubricationModelIter.key()]); - Fs[pair.phase1().index()] += F; - Fs[pair.phase2().index()] -= F; + setF(Fs, pair.phase1().index()) += F; + setF(Fs, pair.phase2().index()) -= F; } return tFs; @@ -581,20 +590,13 @@ Foam::MomentumTransferPhaseSystem::Fs() const template -Foam::autoPtr > -Foam::MomentumTransferPhaseSystem::phiDs +Foam::surfaceScalarField& +Foam::MomentumTransferPhaseSystem::setPhiD ( - const PtrList& rAUs + PtrList& phiDs, const label phasei ) const { - autoPtr > tphiDs - ( - new PtrList(this->phases().size()) - ); - - PtrList& phiDs = tphiDs(); - - forAll(phiDs, phasei) + if (!phiDs.set(phasei)) { phiDs.set ( @@ -603,7 +605,7 @@ Foam::MomentumTransferPhaseSystem::phiDs ( IOobject ( - liftModel::typeName + ":F", + turbulentDispersionModel::typeName + ":phiD", this->mesh_.time().timeName(), this->mesh_, IOobject::NO_READ, @@ -621,6 +623,23 @@ Foam::MomentumTransferPhaseSystem::phiDs ); } + return phiDs[phasei]; +} + + +template +Foam::autoPtr > +Foam::MomentumTransferPhaseSystem::phiDs +( + const PtrList& rAUs +) const +{ + autoPtr > tphiDs + ( + new PtrList(this->phases().size()) + ); + PtrList& phiDs = tphiDs(); + // Add the turbulent dispersion force forAllConstIter ( @@ -638,9 +657,9 @@ Foam::MomentumTransferPhaseSystem::phiDs fvc::snGrad(pair.phase1())*this->mesh_.magSf() ); - phiDs[pair.phase1().index()] += + setPhiD(phiDs, pair.phase1().index()) += fvc::interpolate(rAUs[pair.phase1().index()]*D)*snGradAlpha1; - phiDs[pair.phase2().index()] -= + setPhiD(phiDs, pair.phase2().index()) -= fvc::interpolate(rAUs[pair.phase2().index()]*D)*snGradAlpha1; } diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H index 980f854a9e..2f787ed004 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H @@ -133,6 +133,20 @@ private: //- Turbulent dispersion models turbulentDispersionModelTable turbulentDispersionModels_; + //- Construct element phasei of Fs if not set and return + // Used by Fs() + volVectorField& setF + ( + PtrList& Fs, const label phasei + ) const; + + //- Construct element phasei of phiDs if not set and return + // Used by phiDs() + surfaceScalarField& setPhiD + ( + PtrList& phiDs, const label phasei + ) const; + public: diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H index 136d595418..e23746a3ad 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H @@ -38,22 +38,51 @@ forAll(phases, phasei) PtrList phiFs(phases.size()); { autoPtr > Fs = fluid.Fs(); + + forAll(phases, phasei) + { + phaseModel& phase = phases[phasei]; + + if (Fs().set(phasei)) + { + phiFs.set + ( + phasei, + new surfaceScalarField + ( + IOobject::groupName("phiF", phase.name()), + (fvc::interpolate(rAUs[phasei]*Fs()[phasei]) & mesh.Sf()) + ) + ); + } + } +} +{ autoPtr > phiDs = fluid.phiDs(rAUs); forAll(phases, phasei) { phaseModel& phase = phases[phasei]; - phiFs.set - ( - phasei, - new surfaceScalarField - ( - IOobject::groupName("phiF", phase.name()), - (fvc::interpolate(rAUs[phasei]*Fs()[phasei]) & mesh.Sf()) - + phiDs()[phasei] - ) - ); + if (phiDs().set(phasei)) + { + if (phiFs.set(phasei)) + { + phiFs[phasei] += phiDs()[phasei]; + } + else + { + phiFs.set + ( + phasei, + new surfaceScalarField + ( + IOobject::groupName("phiF", phase.name()), + phiDs()[phasei] + ) + ); + } + } } } @@ -101,12 +130,12 @@ while (pimple.correct()) ghf*fvc::snGrad(rho)*mesh.magSf() ); - PtrList phigs(phases.size()); + PtrList phigFs(phases.size()); forAll(phases, phasei) { phaseModel& phase = phases[phasei]; - phigs.set + phigFs.set ( phasei, ( @@ -118,6 +147,11 @@ while (pimple.correct()) ) ).ptr() ); + + if (phiFs.set(phasei)) + { + phigFs[phasei] += phiFs[phasei]; + } } PtrList phiHbyAs(phases.size()); @@ -177,8 +211,7 @@ while (pimple.correct()) MRF.absolute(phase.phi().oldTime()) - (fvc::interpolate(phase.U().oldTime()) & mesh.Sf()) )/runTime.deltaT() - - phiFs[phasei] - - phigs[phasei] + - phigFs[phasei] ) ); @@ -389,8 +422,7 @@ while (pimple.correct()) + fvc::reconstruct ( alpharAUfs[phasei]*mSfGradp - - phiFs[phasei] - - phigs[phasei] + - phigFs[phasei] ); phase.U().correctBoundaryConditions(); fvOptions.correct(phase.U()); diff --git a/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/constant/phaseProperties b/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/constant/phaseProperties index f7590f98f3..98192e477f 100644 --- a/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/constant/phaseProperties +++ b/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/constant/phaseProperties @@ -151,16 +151,13 @@ heatTransfer ); lift -( -); +(); wallLubrication -( -); +(); turbulentDispersion -( -); +(); // Minimum allowable pressure pMin 10000;