incompressibleDriftFlux: Override the surfaceTensionForce functions

so that the surface tension force is not evaluated or applied as it is
inappropriate for drift-flux modelling of suspensions.

applications/solvers/modules/compressibleMultiphaseVoF/pressureCorrector.C

@@ -64,7 +64,7 @@ void Foam::solvers::compressibleMultiphaseVoF::pressureCorrector()
         const surfaceScalarField phig
         (
             (
-                mixture.surfaceTensionForce(U)
+                surfaceTensionForce()
               - buoyancy.ghf*fvc::snGrad(rho)
             )*rAUf*mesh.magSf()
         );

applications/solvers/modules/compressibleVoF/pressureCorrector.C

@@ -71,7 +71,7 @@ void Foam::solvers::compressibleVoF::pressureCorrector()
         const surfaceScalarField phig
         (
             (
-                interface.surfaceTensionForce()
+                surfaceTensionForce()
               - buoyancy.ghf*fvc::snGrad(rho)
             )*rAUf*mesh.magSf()
         );

applications/solvers/modules/incompressibleDriftFlux/incompressibleDriftFlux.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2022-2023 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2023 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -24,9 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "incompressibleDriftFlux.H"
-#include "localEulerDdtScheme.H"
 #include "fvCorrectPhi.H"
-#include "geometricZeroField.H"
 #include "addToRunTimeSelectionTable.H"
 
 #include "fvmDdt.H"
@@ -45,6 +43,24 @@ namespace solvers
 }
 
 
+// * * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * //
+
+void Foam::solvers::incompressibleDriftFlux::correctInterface()
+{}
+
+
+Foam::tmp<Foam::surfaceScalarField>
+Foam::solvers::incompressibleDriftFlux::surfaceTensionForce() const
+{
+    return surfaceScalarField::New
+    (
+        "surfaceTensionForce",
+        mesh,
+        dimensionedScalar(dimForce/dimVolume, 0)
+    );
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::solvers::incompressibleDriftFlux::incompressibleDriftFlux(fvMesh& mesh)

applications/solvers/modules/incompressibleDriftFlux/incompressibleDriftFlux.H

@@ -147,6 +147,13 @@ protected:
             tmp<volScalarField::Internal>& Sp
         );
 
+        //- Correct the interface properties following mesh-change
+        //  and phase-fraction update
+        virtual void correctInterface();
+
+        //- Return the interface surface tension force for the momentum equation
+        virtual tmp<surfaceScalarField> surfaceTensionForce() const;
+
         //- Return the momentum equation stress term
         virtual tmp<fvVectorMatrix> divDevTau(volVectorField& U)
         {

applications/solvers/modules/incompressibleDriftFlux/pressureCorrector.C

@@ -70,7 +70,7 @@ void Foam::solvers::incompressibleDriftFlux::pressureCorrector()
         surfaceScalarField phig
         (
             (
-                interface.surfaceTensionForce()
+                surfaceTensionForce()
               - buoyancy.ghf*fvc::snGrad(rho)
             )*rAUf*mesh.magSf()
         );

applications/solvers/modules/incompressibleMultiphaseVoF/pressureCorrector.C

@@ -70,7 +70,7 @@ void Foam::solvers::incompressibleMultiphaseVoF::pressureCorrector()
         surfaceScalarField phig
         (
             (
-                mixture.surfaceTensionForce(U)
+                surfaceTensionForce()
               - buoyancy.ghf*fvc::snGrad(rho)
             )*rAUf*mesh.magSf()
         );

applications/solvers/modules/incompressibleVoF/pressureCorrector.C

@@ -70,7 +70,7 @@ void Foam::solvers::incompressibleVoF::pressureCorrector()
         surfaceScalarField phig
         (
             (
-                interface.surfaceTensionForce()
+                surfaceTensionForce()
               - buoyancy.ghf*fvc::snGrad(rho)
             )*rAUf*mesh.magSf()
         );

applications/solvers/modules/twoPhaseVoFSolver/alphaPredictor.C

@@ -227,7 +227,7 @@ void Foam::solvers::twoPhaseVoFSolver::alphaSolve
 
         alpha2 = 1.0 - alpha1;
 
-        interface.correct();
+        correctInterface();
     }
 
     for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
@@ -321,7 +321,7 @@ void Foam::solvers::twoPhaseVoFSolver::alphaSolve
         alpha2 = 1.0 - alpha1;
 
         // Correct only the mixture interface for the interface compression flux
-        interface.correct();
+        correctInterface();
     }
 
     if (alphaApplyPrevCorr && MULESCorr)