From c31789c34cedf1125861ada123d80a67aa1a9b4e Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Wed, 16 Sep 2015 21:29:09 +0100 Subject: [PATCH] reactingEulerFoam: Use PtrListDictionary for list/table of phases This makes looping over the phases much simpler which maintaining support for phase-name lookup. --- .../HeatAndMassTransferPhaseSystem.C | 2 +- .../HeatTransferPhaseSystem.C | 2 +- ...terfaceCompositionPhaseChangePhaseSystem.C | 2 +- .../MomentumTransferPhaseSystem.C | 2 +- .../ThermalPhaseChangePhaseSystem.C | 2 +- .../phaseSystems/phaseSystem/phaseSystem.C | 22 +- .../phaseSystems/phaseSystem/phaseSystem.H | 12 +- .../phaseSystems/phaseSystem/phaseSystemI.H | 4 +- .../phaseSystem/phaseSystemTemplates.H | 2 +- .../reactingMultiphaseEulerFoam/CourantNo.H | 4 +- .../reactingMultiphaseEulerFoam/EEqns.H | 9 +- .../reactingMultiphaseEulerFoam/YEqns.H | 5 +- .../multiphaseSystem/multiphaseSystem.C | 67 ++-- .../reactingMultiphaseEulerFoam/pU/UEqns.H | 8 +- .../reactingMultiphaseEulerFoam/pU/pEqn.H | 77 ++--- .../reactingMultiphaseEulerFoam/pU/pEqnSave.H | 285 ------------------ .../reactingMultiphaseEulerFoam/setRDeltaT.H | 4 +- .../twoPhaseSystem/twoPhaseSystem.C | 4 +- 18 files changed, 93 insertions(+), 420 deletions(-) delete mode 100644 applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqnSave.H diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C index 66d119fd94..bfad310038 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C @@ -270,7 +270,7 @@ Foam::HeatAndMassTransferPhaseSystem::heatTransfer() const forAllConstIter ( - phaseSystem::phaseModelTable, + phaseSystem::phaseModelList, this->phaseModels_, phaseModelIter ) diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.C index 67a3918a32..a0881ee33a 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.C @@ -122,7 +122,7 @@ Foam::HeatTransferPhaseSystem::heatTransfer() const forAllConstIter ( - phaseSystem::phaseModelTable, + phaseSystem::phaseModelList, this->phaseModels_, phaseModelIter ) diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C index 47c994334d..801ceb8caa 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C @@ -71,7 +71,7 @@ massTransfer() const forAllConstIter ( - phaseSystem::phaseModelTable, + phaseSystem::phaseModelList, this->phaseModels_, phaseModelIter ) diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C index cf7fae6f40..0f5470d26e 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C @@ -415,7 +415,7 @@ Foam::MomentumTransferPhaseSystem::momentumTransfer() const forAllConstIter ( - phaseSystem::phaseModelTable, + phaseSystem::phaseModelList, this->phaseModels_, phaseModelIter ) diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C index ffdb41326e..1eb96b3b26 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C @@ -66,7 +66,7 @@ Foam::ThermalPhaseChangePhaseSystem::massTransfer() const forAllConstIter ( - phaseSystem::phaseModelTable, + phaseSystem::phaseModelList, this->phaseModels_, phaseModelIter ) diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C index 31dc2375db..ae23455d88 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C @@ -43,10 +43,10 @@ const Foam::word Foam::phaseSystem::propertiesName("phaseProperties"); Foam::tmp Foam::phaseSystem::calcPhi ( - const phaseModelTable& phaseModels + const phaseModelList& phaseModels ) const { - phaseModelTable::const_iterator phaseModelIter = phaseModels.begin(); + phaseModelList::const_iterator phaseModelIter = phaseModels.begin(); tmp tmpPhi ( @@ -196,7 +196,7 @@ Foam::phaseSystem::~phaseSystem() Foam::tmp Foam::phaseSystem::rho() const { - phaseModelTable::const_iterator phaseModelIter = phaseModels_.begin(); + phaseModelList::const_iterator phaseModelIter = phaseModels_.begin(); tmp tmpRho ( @@ -214,7 +214,7 @@ Foam::tmp Foam::phaseSystem::rho() const Foam::tmp Foam::phaseSystem::U() const { - phaseModelTable::const_iterator phaseModelIter = phaseModels_.begin(); + phaseModelList::const_iterator phaseModelIter = phaseModels_.begin(); tmp tmpU ( @@ -266,7 +266,7 @@ void Foam::phaseSystem::solve() void Foam::phaseSystem::correct() { - forAllIter(phaseModelTable, phaseModels_, phaseModelIter) + forAllIter(phaseModelList, phaseModels_, phaseModelIter) { phaseModelIter().correct(); } @@ -277,7 +277,7 @@ void Foam::phaseSystem::correctKinematics() { bool updateDpdt = false; - forAllIter(phaseModelTable, phaseModels_, phaseModelIter) + forAllIter(phaseModelList, phaseModels_, phaseModelIter) { phaseModelIter().correctKinematics(); @@ -287,14 +287,14 @@ void Foam::phaseSystem::correctKinematics() // Update the pressure time-derivative if required if (updateDpdt) { - dpdt_ = fvc::ddt(phaseModels_.begin()()->thermo().p()); + dpdt_ = fvc::ddt(phaseModels_.begin()().thermo().p()); } } void Foam::phaseSystem::correctThermo() { - forAllIter(phaseModelTable, phaseModels_, phaseModelIter) + forAllIter(phaseModelList, phaseModels_, phaseModelIter) { phaseModelIter().correctThermo(); } @@ -303,7 +303,7 @@ void Foam::phaseSystem::correctThermo() void Foam::phaseSystem::correctTurbulence() { - forAllIter(phaseModelTable, phaseModels_, phaseModelIter) + forAllIter(phaseModelList, phaseModels_, phaseModelIter) { phaseModelIter().correctTurbulence(); } @@ -312,7 +312,7 @@ void Foam::phaseSystem::correctTurbulence() void Foam::phaseSystem::correctEnergyTransport() { - forAllIter(phaseModelTable, phaseModels_, phaseModelIter) + forAllIter(phaseModelList, phaseModels_, phaseModelIter) { phaseModelIter().correctEnergyTransport(); } @@ -325,7 +325,7 @@ bool Foam::phaseSystem::read() { bool readOK = true; - forAllIter(phaseModelTable, phaseModels_, phaseModelIter) + forAllIter(phaseModelList, phaseModels_, phaseModelIter) { readOK &= phaseModelIter().read(); } diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H index bcb268bdc5..95eb7f6a3f 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H @@ -42,7 +42,7 @@ SourceFiles #include "phasePair.H" #include "orderedPhasePair.H" #include "HashPtrTable.H" -#include "PtrDictionary.H" +#include "PtrListDictionary.H" #include "IOMRFZoneList.H" #include "fvIOoptionList.H" @@ -118,7 +118,7 @@ public: > massTransferTable; - typedef PtrDictionary phaseModelTable; + typedef PtrListDictionary phaseModelList; protected: @@ -162,7 +162,7 @@ protected: const fvMesh& mesh_; //- Phase models - phaseModelTable phaseModels_; + phaseModelList phaseModels_; //- Phase pairs phasePairTable phasePairs_; @@ -197,7 +197,7 @@ protected: //- Calculate and return the mixture flux tmp calcPhi ( - const phaseModelTable& phaseModels + const phaseModelList& phaseModels ) const; //- Generate pairs @@ -284,10 +284,10 @@ public: inline const fvMesh& mesh() const; //- Constant access the phase models - inline const phaseModelTable& phases() const; + inline const phaseModelList& phases() const; //- Access the phase models - inline phaseModelTable& phases(); + inline phaseModelList& phases(); //- Constant access the phase pairs inline const phasePairTable& phasePairs() const; diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemI.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemI.H index 343d31f49b..0f9aef5033 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemI.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemI.H @@ -31,14 +31,14 @@ inline const Foam::fvMesh& Foam::phaseSystem::mesh() const } -inline const Foam::phaseSystem::phaseModelTable& +inline const Foam::phaseSystem::phaseModelList& Foam::phaseSystem::phases() const { return phaseModels_; } -inline Foam::phaseSystem::phaseModelTable& +inline Foam::phaseSystem::phaseModelList& Foam::phaseSystem::phases() { return phaseModels_; diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemTemplates.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemTemplates.H index 40d96f9df7..37a3627e1c 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemTemplates.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemTemplates.H @@ -151,7 +151,7 @@ void Foam::phaseSystem::generatePairsAndSubModels HashTable, phasePairKey, phasePairKey::hash> modelTypeTable; - forAllConstIter(phaseModelTable, phaseModels_, phaseModelIter) + forAllConstIter(phaseModelList, phaseModels_, phaseModelIter) { modelTypeTable tempModels; generatePairsAndSubModels diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/CourantNo.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/CourantNo.H index df28acb69f..93c97ac4b0 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/CourantNo.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/CourantNo.H @@ -39,12 +39,12 @@ if (mesh.nInternalFaces()) fvc::surfaceSum(mag(phi))().internalField() ); - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { sumPhi = max ( sumPhi, - fvc::surfaceSum(mag(iter().phi()))().internalField() + fvc::surfaceSum(mag(fluid.phases()[phasei].phi()))().internalField() ); } diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/EEqns.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/EEqns.H index 86150a045d..15a8f6216b 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/EEqns.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/EEqns.H @@ -6,9 +6,10 @@ phaseSystem::heatTransferTable& heatTransfer = heatTransferPtr(); - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; + const volScalarField& alpha = phase; const volScalarField& rho = phase.rho(); const volVectorField& U = phase.U(); @@ -35,9 +36,9 @@ fluid.correctThermo(); -forAllIter(PtrDictionary, fluid.phases(), iter) +forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; Info<< phase.name() << " min/max T " << min(phase.thermo().T()).value() diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/YEqns.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/YEqns.H index 7c4345913d..5a44011fc9 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/YEqns.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/YEqns.H @@ -5,9 +5,10 @@ phaseSystem::massTransferTable& massTransfer(massTransferPtr()); - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; + PtrList& Y = phase.Y(); const volScalarField& alpha = phase; const volScalarField& rho = phase.rho(); diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C index 276dd4e026..cf765d4c07 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C @@ -59,9 +59,9 @@ void Foam::multiphaseSystem::calcAlphas() scalar level = 0.0; alphas_ == 0.0; - forAllIter(PtrDictionary, phases(), iter) + forAll(phases(), i) { - alphas_ += level*iter(); + alphas_ += level*phases()[i]; level += 1.0; } @@ -72,11 +72,9 @@ void Foam::multiphaseSystem::calcAlphas() void Foam::multiphaseSystem::solveAlphas() { PtrList alphaPhiCorrs(phases().size()); - int phasei = 0; - - forAllIter(PtrDictionary, phases(), iter) + forAll(phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = phases()[phasei]; volScalarField& alpha1 = phase; phase.alphaPhi() = @@ -99,9 +97,9 @@ void Foam::multiphaseSystem::solveAlphas() surfaceScalarField& alphaPhiCorr = alphaPhiCorrs[phasei]; - forAllIter(PtrDictionary, phases(), iter2) + forAll(phases(), phasej) { - phaseModel& phase2 = iter2(); + phaseModel& phase2 = phases()[phasej]; volScalarField& alpha2 = phase2; if (&phase2 == &phase) continue; @@ -170,8 +168,6 @@ void Foam::multiphaseSystem::solveAlphas() 0, true ); - - phasei++; } MULES::limitSum(alphaPhiCorrs); @@ -191,10 +187,9 @@ void Foam::multiphaseSystem::solveAlphas() volScalarField divU(fvc::div(fvc::absolute(phi_, phases().first().U()))); - phasei = 0; - forAllIter(PtrDictionary, phases(), iter) + forAll(phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = phases()[phasei]; volScalarField& alpha = phase; surfaceScalarField& alphaPhic = alphaPhiCorrs[phasei]; @@ -245,9 +240,9 @@ void Foam::multiphaseSystem::solveAlphas() } } - forAllConstIter(PtrDictionary, phases(), iter2) + forAll(phases(), phasej) { - const phaseModel& phase2 = iter2(); + const phaseModel& phase2 = phases()[phasej]; const volScalarField& alpha2 = phase2; if (&phase2 == &phase) continue; @@ -294,8 +289,6 @@ void Foam::multiphaseSystem::solveAlphas() << endl; sumAlpha += phase; - - phasei++; } Info<< "Phase-sum volume fraction, min, max = " @@ -506,9 +499,9 @@ Foam::multiphaseSystem::multiphaseSystem 1e-8/pow(average(mesh.V()), 1.0/3.0) ) { - forAllIter(phaseSystem::phaseModelTable, phases(), iter) + forAll(phases(), phasei) { - volScalarField& alphai = iter(); + volScalarField& alphai = phases()[phasei]; mesh.setFluxRequired(alphai.name()); } } @@ -547,9 +540,9 @@ Foam::tmp Foam::multiphaseSystem::surfaceTension ) ); - forAllConstIter(PtrDictionary, phases(), iter) + forAll(phases(), phasej) { - const phaseModel& phase2 = iter(); + const phaseModel& phase2 = phases()[phasej]; if (&phase2 != &phase1) { @@ -591,9 +584,13 @@ Foam::multiphaseSystem::nearInterface() const ) ); - forAllConstIter(PtrDictionary, phases(), iter) + forAll(phases(), phasei) { - tnearInt() = max(tnearInt(), pos(iter() - 0.01)*pos(0.99 - iter())); + tnearInt() = max + ( + tnearInt(), + pos(phases()[phasei] - 0.01)*pos(0.99 - phases()[phasei]) + ); } return tnearInt; @@ -625,10 +622,9 @@ void Foam::multiphaseSystem::solve() PtrList alpha0s(phases().size()); PtrList phiSums(phases().size()); - int phasei = 0; - forAllIter(PtrDictionary, phases(), iter) + forAll(phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = phases()[phasei]; volScalarField& alpha = phase; alpha0s.set @@ -652,8 +648,6 @@ void Foam::multiphaseSystem::solve() dimensionedScalar("0", dimensionSet(0, 3, -1, 0, 0), 0) ) ); - - phasei++; } for @@ -668,18 +662,15 @@ void Foam::multiphaseSystem::solve() { solveAlphas(); - int phasei = 0; - forAllIter(PtrDictionary, phases(), iter) + forAll(phases(), phasei) { - phiSums[phasei] += iter().phi(); - phasei++; + phiSums[phasei] += phases()[phasei].phi(); } } - phasei = 0; - forAllIter(PtrDictionary, phases(), iter) + forAll(phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = phases()[phasei]; volScalarField& alpha = phase; phase.phi() = phiSums[phasei]/nAlphaSubCycles; @@ -691,8 +682,6 @@ void Foam::multiphaseSystem::solve() // Reset the old-time field value alpha.oldTime() = alpha0s[phasei]; alpha.oldTime().timeIndex() = runTime.timeIndex(); - - phasei++; } } else @@ -700,9 +689,9 @@ void Foam::multiphaseSystem::solve() solveAlphas(); } - forAllIter(PtrDictionary, phases(), iter) + forAll(phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = phases()[phasei]; phase.alphaRhoPhi() = fvc::interpolate(phase.rho())*phase.alphaPhi(); } diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/UEqns.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/UEqns.H index 07a5c88194..3b6ee95513 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/UEqns.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/UEqns.H @@ -9,10 +9,10 @@ PtrList UEqns(fluid.phases().size()); phaseSystem::momentumTransferTable& momentumTransfer(momentumTransferPtr()); - int phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; + const volScalarField& alpha = phase; const volScalarField& rho = phase.rho(); volVectorField& U = phase.U(); @@ -31,7 +31,5 @@ PtrList UEqns(fluid.phases().size()); UEqns[phasei].relax(); fvOptions.constrain(UEqns[phasei]); - - phasei++; } } diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H index 442676d1b9..ff759f42c6 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H @@ -2,10 +2,9 @@ PtrList alphafs(fluid.phases().size()); PtrList rAUs(fluid.phases().size()); PtrList alpharAUfs(fluid.phases().size()); -int phasei = 0; -forAllIter(PtrDictionary, fluid.phases(), iter) +forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; const volScalarField& alpha = phase; alphafs.set(phasei, fvc::interpolate(alpha).ptr()); @@ -33,8 +32,6 @@ forAllIter(PtrDictionary, fluid.phases(), iter) fvc::interpolate(max(alpha, phase.residualAlpha())*rAUs[phasei]) ).ptr() ); - - phasei++; } // Turbulent diffusion, particle-pressure, lift and wall-lubrication fluxes @@ -42,10 +39,9 @@ PtrList phiFs(fluid.phases().size()); { autoPtr > Fs = fluid.Fs(); - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; phiFs.set ( @@ -56,8 +52,6 @@ PtrList phiFs(fluid.phases().size()); (fvc::interpolate(rAUs[phasei]*Fs()[phasei]) & mesh.Sf()) ) ); - - phasei++; } } @@ -69,10 +63,9 @@ while (pimple.correct()) PtrList HbyAs(fluid.phases().size()); - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; const volScalarField& alpha = phase; // Correct fixed-flux BCs to be consistent with the velocity BCs @@ -95,8 +88,6 @@ while (pimple.correct()) + max(phase.residualAlpha() - alpha, scalar(0)) *phase.rho()*phase.U().oldTime()/runTime.deltaT() ); - - phasei++; } // Mean density for buoyancy force and p_rgh -> p @@ -109,10 +100,9 @@ while (pimple.correct()) ); PtrList phigs(fluid.phases().size()); - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; phigs.set ( @@ -126,8 +116,6 @@ while (pimple.correct()) ) ).ptr() ); - - phasei++; } PtrList phiHbyAs(fluid.phases().size()); @@ -146,10 +134,9 @@ while (pimple.correct()) dimensionedScalar("phiHbyA", dimArea*dimVelocity, 0) ); - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; const volScalarField& alpha = phase; // ddtPhiCorr filter -- only apply in pure(ish) phases @@ -223,8 +210,6 @@ while (pimple.correct()) } phiHbyA += alphafs[phasei]*phiHbyAs[phasei]; - - phasei++; } MRF.makeRelative(phiHbyA); @@ -242,11 +227,9 @@ while (pimple.correct()) dimensionedScalar("rAUf", dimensionSet(-1, 3, 1, 0, 0), 0) ); - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { rAUf += alphafs[phasei]*alpharAUfs[phasei]; - phasei++; } rAUf = mag(rAUf); @@ -255,12 +238,10 @@ while (pimple.correct()) { surfaceScalarField::GeometricBoundaryField phib(phi.boundaryField()); phib = 0; - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; phib += alphafs[phasei].boundaryField()*phase.phi().boundaryField(); - phasei++; } setSnGrad @@ -273,10 +254,9 @@ while (pimple.correct()) } PtrList pEqnComps(fluid.phases().size()); - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; if (phase.compressible()) { @@ -338,8 +318,6 @@ while (pimple.correct()) ); } } - - phasei++; } // Cache p prior to solve for density update @@ -358,17 +336,14 @@ while (pimple.correct()) { fvScalarMatrix pEqn(pEqnIncomp); - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; if (phase.compressible()) { pEqn += pEqnComps[phasei]; } - - phasei++; } solve @@ -385,10 +360,9 @@ while (pimple.correct()) surfaceScalarField mSfGradp("mSfGradp", pEqnIncomp.flux()/rAUf); - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; phase.phi() = phiHbyAs[phasei] + alpharAUfs[phasei]*mSfGradp; @@ -397,8 +371,6 @@ while (pimple.correct()) { phase.divU(-pEqnComps[phasei] & p_rgh); } - - phasei++; } // Optionally relax pressure for velocity correction @@ -406,10 +378,9 @@ while (pimple.correct()) mSfGradp = pEqnIncomp.flux()/rAUf; - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; phase.U() = HbyAs[phasei] @@ -421,8 +392,6 @@ while (pimple.correct()) ); phase.U().correctBoundaryConditions(); fvOptions.correct(phase.U()); - - phasei++; } } } @@ -434,9 +403,9 @@ while (pimple.correct()) p_rgh = p - rho*gh; // Update densities from change in p_rgh - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - phaseModel& phase = iter(); + phaseModel& phase = fluid.phases()[phasei]; phase.rho()() += phase.thermo().psi()*(p_rgh - p_rgh_0); } diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqnSave.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqnSave.H deleted file mode 100644 index e94f68d8f2..0000000000 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqnSave.H +++ /dev/null @@ -1,285 +0,0 @@ -// --- Pressure corrector loop -while (pimple.correct()) -{ - // rho1 = rho10 + psi1*p_rgh; - // rho2 = rho20 + psi2*p_rgh; - - // tmp pEqnComp1; - // tmp pEqnComp2; - - // //if (transonic) - // //{ - // //} - // //else - // { - // surfaceScalarField phid1("phid1", fvc::interpolate(psi1)*phi1); - // surfaceScalarField phid2("phid2", fvc::interpolate(psi2)*phi2); - - // pEqnComp1 = - // fvc::ddt(rho1) + psi1*correction(fvm::ddt(p_rgh)) - // + fvc::div(phid1, p_rgh) - // - fvc::Sp(fvc::div(phid1), p_rgh); - - // pEqnComp2 = - // fvc::ddt(rho2) + psi2*correction(fvm::ddt(p_rgh)) - // + fvc::div(phid2, p_rgh) - // - fvc::Sp(fvc::div(phid2), p_rgh); - // } - - PtrList alphafs(fluid.phases().size()); - PtrList HbyAs(fluid.phases().size()); - PtrList phiHbyAs(fluid.phases().size()); - PtrList rAUs(fluid.phases().size()); - PtrList rAlphaAUfs(fluid.phases().size()); - - int phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) - { - phaseModel& phase = iter(); - - MRF.makeAbsolute(phase.phi().oldTime()); - MRF.makeAbsolute(phase.phi()); - - HbyAs.set(phasei, new volVectorField(phase.U())); - phiHbyAs.set(phasei, new surfaceScalarField(1.0*phase.phi())); - - phasei++; - } - - surfaceScalarField phiHbyA - ( - IOobject - ( - "phiHbyA", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - mesh, - dimensionedScalar("phiHbyA", dimArea*dimVelocity, 0) - ); - - volScalarField rho("rho", fluid.rho()); - surfaceScalarField ghSnGradRho(ghf*fvc::snGrad(rho)*mesh.magSf()); - - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) - { - phaseModel& phase = iter(); - const volScalarField& alpha = phase; - - alphafs.set(phasei, fvc::interpolate(alpha).ptr()); - alphafs[phasei].rename("hmm" + alpha.name()); - - volScalarField dragCoeffi - ( - IOobject - ( - "dragCoeffi", - runTime.timeName(), - mesh - ), - fluid.Kd(phase), - zeroGradientFvPatchScalarField::typeName - ); - dragCoeffi.correctBoundaryConditions(); - - rAUs.set(phasei, (1.0/(UEqns[phasei].A() + dragCoeffi)).ptr()); - rAlphaAUfs.set - ( - phasei, - ( - alphafs[phasei]*fvc::interpolate(phase.rho()) - /fvc::interpolate(UEqns[phasei].A() + dragCoeffi) - ).ptr() - ); - - HbyAs[phasei] = rAUs[phasei]*UEqns[phasei].H(); - - phiHbyAs[phasei] = - ( - (fvc::interpolate(HbyAs[phasei]) & mesh.Sf()) - + rAlphaAUfs[phasei]*fvc::ddtCorr(phase.U(), phase.phi()) - ); - MRF.makeRelative(phiHbyAs[phasei]); - MRF.makeRelative(phase.phi().oldTime()); - MRF.makeRelative(phase.phi()); - - phiHbyAs[phasei] += - rAlphaAUfs[phasei] - *( - fluid.surfaceTension(phase)*mesh.magSf() - + (fvc::interpolate(phase.rho() - rho))*(g & mesh.Sf()) - - ghSnGradRho - )/fvc::interpolate(phase.rho()); - - forAllConstIter - ( - phaseSystem::KdTable, - fluid.Kds(), - KdIter - ) - { - const volScalarField& K(*KdIter()); - - const phasePair& pair(fluid.phasePairs()[KdIter.key()]); - - const phaseModel* phase1 = &pair.phase1(); - const phaseModel* phase2 = &pair.phase2(); - - forAllConstIter(phasePair, pair, iter) - { - if (phase1 == &phase) - { - phiHbyAs[phasei] += - fvc::interpolate(K) - /fvc::interpolate(UEqns[phasei].A() + dragCoeffi) - *phase2->phi(); - - HbyAs[phasei] += rAUs[phasei]*K*phase2->U(); - } - - Swap(phase1, phase2); - } - } - - phiHbyA += alphafs[phasei]*phiHbyAs[phasei]; - - phasei++; - } - - surfaceScalarField rAUf - ( - IOobject - ( - "rAUf", - runTime.timeName(), - mesh - ), - mesh, - dimensionedScalar("rAUf", dimensionSet(-1, 3, 1, 0, 0), 0) - ); - - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) - { - phaseModel& phase = iter(); - rAUf += mag(alphafs[phasei]*rAlphaAUfs[phasei]) - /fvc::interpolate(phase.rho()); - - phasei++; - } - - // Update the fixedFluxPressure BCs to ensure flux consistency - { - surfaceScalarField::GeometricBoundaryField phib(phi.boundaryField()); - phib = 0; - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) - { - phaseModel& phase = iter(); - - phib += - alphafs[phasei].boundaryField() - *(mesh.Sf().boundaryField() & phase.U().boundaryField()); - - phasei++; - } - - setSnGrad - ( - p_rgh.boundaryField(), - ( - phiHbyA.boundaryField() - MRF.relative(phib) - )/(mesh.magSf().boundaryField()*rAUf.boundaryField()) - ); - } - - while (pimple.correctNonOrthogonal()) - { - fvScalarMatrix pEqnIncomp - ( - fvc::div(phiHbyA) - - fvm::laplacian(rAUf, p_rgh) - ); - - pEqnIncomp.setReference(pRefCell, pRefValue); - - solve - ( - // ( - // (alpha1/rho1)*pEqnComp1() - // + (alpha2/rho2)*pEqnComp2() - // ) + - pEqnIncomp, - mesh.solver(p_rgh.select(pimple.finalInnerIter())) - ); - - if (pimple.finalNonOrthogonalIter()) - { - surfaceScalarField mSfGradp("mSfGradp", pEqnIncomp.flux()/rAUf); - - phasei = 0; - phi = dimensionedScalar("phi", phi.dimensions(), 0); - forAllIter(PtrDictionary, fluid.phases(), iter) - { - phaseModel& phase = iter(); - - phase.phi() = - phiHbyAs[phasei] - + rAlphaAUfs[phasei]*mSfGradp/fvc::interpolate(phase.rho()); - phi += - alphafs[phasei]*phiHbyAs[phasei] - + mag(alphafs[phasei]*rAlphaAUfs[phasei]) - *mSfGradp/fvc::interpolate(phase.rho()); - - phasei++; - } - - // dgdt = - - // ( - // pos(alpha2)*(pEqnComp2 & p)/rho2 - // - pos(alpha1)*(pEqnComp1 & p)/rho1 - // ); - - p_rgh.relax(); - - p = p_rgh + rho*gh; - - mSfGradp = pEqnIncomp.flux()/rAUf; - - phasei = 0; - forAllIter(PtrDictionary, fluid.phases(), iter) - { - phaseModel& phase = iter(); - - phase.U() = - HbyAs[phasei] - + fvc::reconstruct - ( - rAlphaAUfs[phasei] - *( - fvc::interpolate(phase.rho() - rho) - *(g & mesh.Sf()) - - ghSnGradRho - + mSfGradp - ) - )/phase.rho(); - - phase.U().correctBoundaryConditions(); - - phasei++; - } - } - } - - //p = max(p, pMin); - - // rho1 = rho10 + psi1*p_rgh; - // rho2 = rho20 + psi2*p_rgh; - - // Dp1Dt = fvc::DDt(phi1, p_rgh); - // Dp2Dt = fvc::DDt(phi2, p_rgh); -} diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/setRDeltaT.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/setRDeltaT.H index 3a609e41a8..aae1d9072c 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/setRDeltaT.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/setRDeltaT.H @@ -12,9 +12,9 @@ surfaceScalarField maxPhi("maxPhi", phi); - forAllIter(PtrDictionary, fluid.phases(), iter) + forAll(fluid.phases(), phasei) { - maxPhi = max(maxPhi, mag(iter().phi())); + maxPhi = max(maxPhi, mag(fluid.phases()[phasei].phi())); } // Set the reciprocal time-step from the local Courant number diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C index 13a98b885f..401c4f622c 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C @@ -57,8 +57,8 @@ Foam::twoPhaseSystem::twoPhaseSystem ) : phaseSystem(mesh), - phase1_(phaseModels_.first()), - phase2_(phaseModels_.last()) + phase1_(phaseModels_[0]), + phase2_(phaseModels_[1]) { phase2_.volScalarField::operator=(scalar(1) - phase1_);