diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C index 9646a6a80c1c167f9eebf08285ee963778224507..a209933b0950b8bc9572737b01167fbca22cdfd0 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C @@ -155,6 +155,63 @@ Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::Kdf } +template<class BasePhaseSystem> +Foam::tmp<Foam::volScalarField> +Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::Kd +( + const Foam::phaseModel& phase +) const +{ + tmp<volScalarField> tKd + ( + new volScalarField + ( + IOobject + ( + IOobject::groupName("Kd", phase.name()), + this->mesh_.time().timeName(), + this->mesh_ + ), + this->mesh_, + dimensionedScalar + ( + IOobject::groupName("Kd", phase.name()), + dimensionSet(1, -3, -1, 0, 0), + 0 + ) + ) + ); + + // Add the implicit part of the drag force + forAllConstIter + ( + phaseSystem::KdTable, + Kds_, + KdIter + ) + { + const volScalarField& K(*KdIter()); + + const phasePair& pair(this->phasePairs_[KdIter.key()]); + + const phaseModel* phase1 = &pair.phase1(); + const phaseModel* phase2 = &pair.phase2(); + + forAllConstIter(phasePair, pair, iter) + { + if (phase1 == &phase) + { + tKd() += K; + } + + Swap(phase1, phase2); + } + } + + return tKd; +} + + template<class BasePhaseSystem> Foam::tmp<Foam::volScalarField> Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::Vm @@ -385,6 +442,7 @@ Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::momentumTransfer() const } // Add the implicit part of the drag force + /* ***HGW Currently this is handled in the pEqn forAllConstIter ( phaseSystem::KdTable, @@ -401,13 +459,14 @@ Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::momentumTransfer() const forAllConstIter(phasePair, pair, iter) { - const volVectorField& U(phase->U()); + const volVectorField& U = phase->U(); *eqns[phase->name()] -= fvm::Sp(K, U); Swap(phase, otherPhase); } } + */ // Update the virtual mass coefficients forAllConstIter @@ -437,8 +496,8 @@ Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::momentumTransfer() const forAllConstIter(phasePair, pair, iter) { - const volVectorField& U(phase->U()); - const surfaceScalarField& phi(phase->phi()); + const volVectorField& U = phase->U(); + const surfaceScalarField& phi = phase->phi(); *eqns[phase->name()] += - Vm diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H index 5c9f057bedc2fa2f8871325a735e4416cad4ce00..b6aeadd34ca656cf7b7748cf52d01d73148519bc 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H @@ -148,12 +148,21 @@ public: // Member Functions + //- Constant access to drag coefficients + virtual const phaseSystem::KdTable& Kds() const + { + return Kds_; + } + //- Return the drag coefficient virtual tmp<volScalarField> Kd(const phasePairKey& key) const; //- Return the face drag coefficient virtual tmp<surfaceScalarField> Kdf(const phasePairKey& key) const; + //- Return the drag coefficient for phase + virtual tmp<volScalarField> Kd(const phaseModel& phase) const; + //- Return the virtual mass coefficient virtual tmp<volScalarField> Vm(const phasePairKey& key) const; diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.C index d64a5c862ef04bc5fba975dd38a9ba3eecaec66b..090d29fafe00eecf3704bc714eacbd4f7c31aecc 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/MovingPhaseModel/MovingPhaseModel.C @@ -260,6 +260,7 @@ Foam::MovingPhaseModel<BasePhaseModel>::UEqn() fvm::ddt(*this, this->thermo().rho(), U_) + fvm::div(alphaRhoPhi_, U_) - fvm::Sp(continuityError_, U_) + + this->fluid().MRF().DDt(*this*this->thermo().rho(), U_) + turbulence_->divDevRhoReff(U_) ); } diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C index 3c5dafdfa6d8f4a08c51e95d0ea883484f69730a..dee40846b63f6a67bb9882ffa0bc5afd0ca7f9df 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.C @@ -282,8 +282,6 @@ void Foam::phaseSystem::correctKinematics() updateDpdt = updateDpdt || phaseModelIter().thermo().dpdt(); } - //phaseModelTable::iterator iter = phaseModels_.begin(); - // Update the pressure time-derivative if required if (updateDpdt) { diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H index 84c69aa655c64e2583d082968e12a18d35ae4a5f..bcb268bdc5f37df76c4a19d2df38226e9549e553 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystem.H @@ -118,10 +118,6 @@ public: > massTransferTable; - // typedef - // HashTable<autoPtr<phaseModel>, word, word::hash> - // phaseModelTable; - typedef PtrDictionary<phaseModel> phaseModelTable; @@ -284,48 +280,56 @@ public: // Member Functions + //- Constant access the mesh + inline const fvMesh& mesh() const; + + //- Constant access the phase models + inline const phaseModelTable& phases() const; + + //- Access the phase models + inline phaseModelTable& phases(); + + //- Constant access the phase pairs + inline const phasePairTable& phasePairs() const; + //- Return the mixture density tmp<volScalarField> rho() const; //- Return the mixture velocity tmp<volVectorField> U() const; - //- Return the surface tension coefficient - tmp<volScalarField> sigma(const phasePairKey& key) const; - - //- Return the drag coefficient - virtual tmp<volScalarField> Kd(const phasePairKey& key) const = 0; - - //- Return the face drag coefficient - virtual tmp<surfaceScalarField> Kdf(const phasePairKey& key) const = 0; + //- Constant access the mixture flux + inline const surfaceScalarField& phi() const; - //- Return the virtual mass coefficient - virtual tmp<volScalarField> Vm(const phasePairKey& key) const = 0; + //- Access the mixture flux + inline surfaceScalarField& phi(); - //- Return the face virtual mass coefficient - virtual tmp<surfaceScalarField> Vmf(const phasePairKey& key) const = 0; + //- Constant access the rate of change of the pressure + inline const volScalarField& dpdt() const; - //- Return the combined force (lift + wall-lubrication) - virtual tmp<volVectorField> F(const phasePairKey& key) const = 0; + //- Access the rate of change of the pressure + inline volScalarField& dpdt(); - //- Return the combined face-force (lift + wall-lubrication) - virtual tmp<surfaceScalarField> Ff(const phasePairKey& key) const = 0; - - //- Return the turbulent diffusivity - // Multiplies the phase-fraction gradient - virtual tmp<volScalarField> D(const phasePairKey& key) const = 0; + //- Return the surface tension coefficient + tmp<volScalarField> sigma(const phasePairKey& key) const; - //- Return the interfacial mass flow rate - virtual tmp<volScalarField> dmdt(const phasePairKey& key) const = 0; + //- Return MRF zones + inline const IOMRFZoneList& MRF() const; - //- Return the momentum transfer matrices - virtual autoPtr<momentumTransferTable> momentumTransfer() const = 0; + //- Optional FV-options + inline fv::IOoptionList& fvOptions() const; - //- Return the heat transfer matrices - virtual autoPtr<heatTransferTable> heatTransfer() const = 0; + //- Access a sub model between a phase pair + template <class modelType> + const modelType& lookupSubModel(const phasePair& key) const; - //- Return the mass transfer matrices - virtual autoPtr<massTransferTable> massTransfer() const = 0; + //- Access a sub model between two phases + template <class modelType> + const modelType& lookupSubModel + ( + const phaseModel& dispersed, + const phaseModel& continuous + ) const; //- Solve for the phase fractions virtual void solve(); @@ -347,48 +351,6 @@ public: //- Read base phaseProperties dictionary virtual bool read(); - - - // Access - - //- Constant access the mesh - inline const fvMesh& mesh() const; - - //- Constant access the phase models - inline const phaseModelTable& phases() const; - - //- Access the phase models - inline phaseModelTable& phases(); - - //- Constant access the mixture flux - inline const surfaceScalarField& phi() const; - - //- Access the mixture flux - inline surfaceScalarField& phi(); - - //- Constant access the rate of change of the pressure - inline const volScalarField& dpdt() const; - - //- Access the rate of change of the pressure - inline volScalarField& dpdt(); - - //- Return MRF zones - inline const IOMRFZoneList& MRF() const; - - //- Optional FV-options - inline fv::IOoptionList& fvOptions() const; - - //- Access a sub model between a phase pair - template <class modelType> - const modelType& lookupSubModel(const phasePair& key) const; - - //- Access a sub model between two phases - template <class modelType> - const modelType& lookupSubModel - ( - const phaseModel& dispersed, - const phaseModel& continuous - ) const; }; diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemI.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemI.H index 35b61ab4516a4d4695ad39ed0bb4ad81aa19d8c3..343d31f49b7d649dcb1d7142362b0a1614e98ebe 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemI.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseSystem/phaseSystemI.H @@ -45,6 +45,13 @@ Foam::phaseSystem::phases() } +inline const Foam::phaseSystem::phasePairTable& +Foam::phaseSystem::phasePairs() const +{ + return phasePairs_; +} + + inline const Foam::surfaceScalarField& Foam::phaseSystem::phi() const { return phi_; diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/UEqns.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/UEqns.H index 31bd20ff7a6e8a8f6aeda844aebed61d74f1f9df..9252630429081e7755931ebf9f4ae466f7506a09 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/UEqns.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pU/UEqns.H @@ -1,8 +1,5 @@ Info<< "Constructing momentum equations" << endl; -MRF.correctBoundaryVelocity(U1); -MRF.correctBoundaryVelocity(U2); - fvVectorMatrix U1Eqn(U1, rho1.dimensions()*U1.dimensions()*dimVol/dimTime); fvVectorMatrix U2Eqn(U2, rho2.dimensions()*U2.dimensions()*dimVol/dimTime); diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/UEqns.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/UEqns.H index eb1456bc15807fc77b930c42de58ea1df68b6aa4..83c3a9ce7f736029755385c297c4834c9dc205a0 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/UEqns.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/pUf/UEqns.H @@ -1,8 +1,5 @@ Info<< "Constructing face momentum equations" << endl; -MRF.correctBoundaryVelocity(U1); -MRF.correctBoundaryVelocity(U2); - fvVectorMatrix U1Eqn(U1, rho1.dimensions()*U1.dimensions()*dimVol/dimTime); fvVectorMatrix U2Eqn(U2, rho2.dimensions()*U2.dimensions()*dimVol/dimTime); diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C index 50b79b5672ba642d11b0a867d1698b147ca144ad..13a98b885f9b0af9730590f543ec303d6a497c8e 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.C @@ -85,6 +85,12 @@ Foam::twoPhaseSystem::sigma() const } +const Foam::dragModel& Foam::twoPhaseSystem::drag(const phaseModel& phase) const +{ + return lookupSubModel<dragModel>(phase, otherPhase(phase)); +} + + Foam::tmp<Foam::volScalarField> Foam::twoPhaseSystem::Kd() const { @@ -105,6 +111,13 @@ Foam::twoPhaseSystem::Kdf() const } +const Foam::virtualMassModel& +Foam::twoPhaseSystem::virtualMass(const phaseModel& phase) const +{ + return lookupSubModel<virtualMassModel>(phase, otherPhase(phase)); +} + + Foam::tmp<Foam::volScalarField> Foam::twoPhaseSystem::Vm() const { @@ -409,17 +422,4 @@ void Foam::twoPhaseSystem::solve() } -const Foam::dragModel& Foam::twoPhaseSystem::drag(const phaseModel& phase) const -{ - return lookupSubModel<dragModel>(phase, otherPhase(phase)); -} - - -const Foam::virtualMassModel& -Foam::twoPhaseSystem::virtualMass(const phaseModel& phase) const -{ - return lookupSubModel<virtualMassModel>(phase, otherPhase(phase)); -} - - // ************************************************************************* // diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.H index 0ac971b790e328735fa56fec79a25276facb52d5..deca55fa6e530e72d3de74c72f315c59f6339d2a 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.H +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.H @@ -53,6 +53,35 @@ class twoPhaseSystem : public phaseSystem { + // Private member functions + + //- Return the drag coefficient for phase pair + virtual tmp<volScalarField> Kd(const phasePairKey& key) const = 0; + + //- Return the face drag coefficient for phase pair + virtual tmp<surfaceScalarField> Kdf(const phasePairKey& key) const = 0; + + //- Return the virtual mass coefficient for phase pair + virtual tmp<volScalarField> Vm(const phasePairKey& key) const = 0; + + //- Return the face virtual mass coefficient for phase pair + virtual tmp<surfaceScalarField> Vmf(const phasePairKey& key) const = 0; + + //- Return the combined force (lift + wall-lubrication) for phase pair + virtual tmp<volVectorField> F(const phasePairKey& key) const = 0; + + //- Return the combined face-force (lift + wall-lubrication) + // for phase pair + virtual tmp<surfaceScalarField> Ff(const phasePairKey& key) const = 0; + + //- Return the turbulent diffusivity for phase pair + // Multiplies the phase-fraction gradient + virtual tmp<volScalarField> D(const phasePairKey& key) const = 0; + + //- Return the interfacial mass flow rate for phase pair + virtual tmp<volScalarField> dmdt(const phasePairKey& key) const = 0; + + protected: // Protected data @@ -103,28 +132,50 @@ public: // Member Functions - //- Solve for the phase fractions - virtual void solve(); + //- Constant access phase model 1 + const phaseModel& phase1() const; + + //- Access phase model 1 + phaseModel& phase1(); + + //- Constant access phase model 2 + const phaseModel& phase2() const; + + //- Access phase model 2 + phaseModel& phase2(); + + //- Constant access the phase not given as an argument + const phaseModel& otherPhase + ( + const phaseModel& phase + ) const; + + //- Return the momentum transfer matrices + virtual autoPtr<momentumTransferTable> momentumTransfer() const = 0; + + //- Return the heat transfer matrices + virtual autoPtr<heatTransferTable> heatTransfer() const = 0; + + //- Return the mass transfer matrices + virtual autoPtr<massTransferTable> massTransfer() const = 0; using phaseSystem::sigma; - using phaseSystem::Kd; - using phaseSystem::Kdf; - using phaseSystem::Vm; - using phaseSystem::Vmf; - using phaseSystem::F; - using phaseSystem::Ff; - using phaseSystem::D; - using phaseSystem::dmdt; //- Return the surface tension coefficient tmp<volScalarField> sigma() const; + //- Return the drag model for the given phase + const dragModel& drag(const phaseModel& phase) const; + //- Return the drag coefficient tmp<volScalarField> Kd() const; //- Return the face drag coefficient tmp<surfaceScalarField> Kdf() const; + //- Return the virtual mass model for the given phase + const virtualMassModel& virtualMass(const phaseModel& phase) const; + //- Return the virtual mass coefficient tmp<volScalarField> Vm() const; @@ -144,32 +195,8 @@ public: //- Return the interfacial mass flow rate tmp<volScalarField> dmdt() const; - - // Access - - //- Constant access phase model 1 - const phaseModel& phase1() const; - - //- Access phase model 1 - phaseModel& phase1(); - - //- Constant access phase model 2 - const phaseModel& phase2() const; - - //- Access phase model 2 - phaseModel& phase2(); - - //- Constant access the phase not given as an argument - const phaseModel& otherPhase - ( - const phaseModel& phase - ) const; - - //- Return the drag model for the given phase - const dragModel& drag(const phaseModel& phase) const; - - //- Return the virtual mass model for the given phase - const virtualMassModel& virtualMass(const phaseModel& phase) const; + //- Solve for the phase fractions + virtual void solve(); };