reactingMultiphaseEulerFoam: Added support for turbulent dispersion

OpenFOAM · Sep 17, 2015 · a1f2f592989cb46bf56946acc501860abf9d1051 · a1f2f59
1 parent 9f666b0
commit a1f2f592989cb46bf56946acc501860abf9d1051
diff --git a/...acialModels/turbulentDispersionModels/turbulentDispersionModel/turbulentDispersionModel.C b/...acialModels/turbulentDispersionModels/turbulentDispersionModel/turbulentDispersionModel.C
@@ -26,6 +26,8 @@ License
 #include "turbulentDispersionModel.H"
 #include "phasePair.H"
 #include "fvcGrad.H"
+#include "surfaceInterpolate.H"
+#include "fvcSnGrad.H"

 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

@@ -66,4 +68,11 @@ Foam::turbulentDispersionModel::F() const
 }


+Foam::tmp<Foam::surfaceScalarField>
+Foam::turbulentDispersionModel::Ff() const
+{
+    return fvc::interpolate(D())*fvc::snGrad(pair_.dispersed());
+}
+
+
 // ************************************************************************* //
diff --git a/...acialModels/turbulentDispersionModels/turbulentDispersionModel/turbulentDispersionModel.H b/...acialModels/turbulentDispersionModels/turbulentDispersionModel/turbulentDispersionModel.H
@@ -120,6 +120,9 @@ public:

        //- Turbulent dispersion force
        virtual tmp<volVectorField> F() const;
+
+        //- Turbulent dispersion force on faces
+        virtual tmp<surfaceScalarField> Ff() const;
 };


diff --git a/...rFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C b/...rFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C
@@ -38,6 +38,7 @@ License
 #include "fvmDiv.H"
 #include "fvmSup.H"
 #include "fvcDiv.H"
+#include "fvcSnGrad.H"
 #include "fvMatrix.H"

 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
@@ -575,24 +576,75 @@ Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::Fs() const
        Fs[pair.phase2().index()] -= F;
    }

+    return tFs;
+}
+
+
+template<class BasePhaseSystem>
+Foam::autoPtr<Foam::PtrList<Foam::surfaceScalarField> >
+Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::phiDs
+(
+    const PtrList<volScalarField>& rAUs
+) const
+{
+    autoPtr<PtrList<surfaceScalarField> > tphiDs
+    (
+        new PtrList<surfaceScalarField>(this->phases().size())
+    );
+
+    PtrList<surfaceScalarField>& phiDs = tphiDs();
+
+    forAll(phiDs, phasei)
+    {
+        phiDs.set
+        (
+            phasei,
+            new surfaceScalarField
+            (
+                IOobject
+                (
+                    liftModel::typeName + ":F",
+                    this->mesh_.time().timeName(),
+                    this->mesh_,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE,
+                    false
+                ),
+                this->mesh_,
+                dimensionedScalar
+                (
+                    "zero",
+                    dimTime*dimArea*turbulentDispersionModel::dimF/dimDensity,
+                    0
+                )
+            )
+        );
+    }
+
    // Add the turbulent dispersion force
-    // forAllConstIter
-    // (
-    //     turbulentDispersionModelTable,
-    //     turbulentDispersionModels_,
-    //     turbulentDispersionModelIter
-    // )
-    // {
-    //     const volVectorField F(turbulentDispersionModelIter()->F<vector>());
-
-    //     const phasePair&
-    //         pair(this->phasePairs_[turbulentDispersionModelIter.key()]);
-
-    //     *eqns[pair.phase1().name()] += F;
-    //     *eqns[pair.phase2().name()] -= F;
-    // }
+    forAllConstIter
+    (
+        turbulentDispersionModelTable,
+        turbulentDispersionModels_,
+        turbulentDispersionModelIter
+    )
+    {
+        const phasePair&
+            pair(this->phasePairs_[turbulentDispersionModelIter.key()]);

-    return tFs;
+        const volScalarField D(turbulentDispersionModelIter()->D());
+        const surfaceScalarField snGradAlpha1
+        (
+            fvc::snGrad(pair.phase1())*this->mesh_.magSf()
+        );
+
+        phiDs[pair.phase1().index()] +=
+            fvc::interpolate(rAUs[pair.phase1().index()]*D)*snGradAlpha1;
+        phiDs[pair.phase2().index()] -=
+            fvc::interpolate(rAUs[pair.phase2().index()]*D)*snGradAlpha1;
+    }
+
+    return tphiDs;
 }


diff --git a/...rFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H b/...rFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.H
@@ -175,6 +175,12 @@ public:
        //- Return the combined force (lift + wall-lubrication)
        virtual autoPtr<PtrList<volVectorField> > Fs() const;

+        //- Return the turbulent dispersion force on faces for phase pair
+        virtual autoPtr<PtrList<surfaceScalarField> > phiDs
+        (
+            const PtrList<volScalarField>& rAUs
+        ) const;
+
        //- Return the combined face-force (lift + wall-lubrication)
        virtual tmp<surfaceScalarField> Ff(const phasePairKey& key) const;

diff --git a/...tiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.H b/...tiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystem.H
@@ -169,6 +169,12 @@ public:
        //- Return the combined force (lift + wall-lubrication) for phase pair
        virtual autoPtr<PtrList<Foam::volVectorField> > Fs() const = 0;

+        //- Return the turbulent dispersion force on faces for phase pair
+        virtual autoPtr<PtrList<Foam::surfaceScalarField> > phiDs
+        (
+            const PtrList<volScalarField>& rAUs
+        ) const = 0;
+
        //- Return the total interfacial mass transfer rate for phase
        virtual tmp<volScalarField> dmdt(const phaseModel& phase) const = 0;

diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H b/applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/pU/pEqn.H
@@ -34,10 +34,11 @@ forAll(phases, phasei)
    );
 }

-// Turbulent diffusion, particle-pressure, lift and wall-lubrication fluxes
+// Lift, wall-lubrication and turbulent diffusion fluxes
 PtrList<surfaceScalarField> phiFs(phases.size());
 {
    autoPtr<PtrList<volVectorField> > Fs = fluid.Fs();
+    autoPtr<PtrList<surfaceScalarField> > phiDs = fluid.phiDs(rAUs);

    forAll(phases, phasei)
    {
@@ -50,6 +51,7 @@ PtrList<surfaceScalarField> phiFs(phases.size());
            (
                IOobject::groupName("phiF", phase.name()),
                (fvc::interpolate(rAUs[phasei]*Fs()[phasei]) & mesh.Sf())
+              + phiDs()[phasei]
            )
        );
    }