cavitationModels: Made insensitive to phase order

The cavitation models used by the interFoam solver and the
compressibleVoF solver module can now be applied regardless of the
ordering of the liquid and vapour phases. A "liquid" keyword is now
required in the model specification in order to control which phase is
considered to be the condensed liquid state. Previously the liquid phase
was assumed to be the first of the two phases.

applications/solvers/modules/fluid/compressibleVoF/fvModels/VoFCavitation/VoFCavitation.C

@@ -96,21 +96,21 @@ void Foam::fv::compressible::VoFCavitation::addSup
 
     if (fieldName == alphaName_)
     {
-        volScalarField::Internal alphalCoeff
+        const volScalarField::Internal alpha1Coeff
         (
             1.0/mixture_.rho1()()
           - mixture_.alpha1()()*(1.0/mixture_.rho1()() - 1.0/mixture_.rho2()())
         );
 
-        const Pair<tmp<volScalarField::Internal>> mDotAlphal
+        const Pair<tmp<volScalarField::Internal>> mDot12Alpha
         (
-            cavitation_->mDotAlphal()
+            cavitation_->mDot12Alpha()
         );
 
-        const volScalarField::Internal vDotcAlphal(alphalCoeff*mDotAlphal[0]());
+        const volScalarField::Internal vDot1Alpha(alpha1Coeff*mDot12Alpha[0]());
-        const volScalarField::Internal vDotvAlphal(alphalCoeff*mDotAlphal[1]());
+        const volScalarField::Internal vDot2Alpha(alpha1Coeff*mDot12Alpha[1]());
 
-        eqn += fvm::Sp(vDotvAlphal - vDotcAlphal, eqn.psi()) + vDotcAlphal;
+        eqn += fvm::Sp(-vDot2Alpha - vDot1Alpha, eqn.psi()) + vDot1Alpha;
     }
 }
 
@@ -140,17 +140,17 @@ void Foam::fv::compressible::VoFCavitation::addSup
             1.0/mixture_.rho1()() - 1.0/mixture_.rho2()()
         );
 
-        const Pair<tmp<volScalarField::Internal>> mDotP
+        const Pair<tmp<volScalarField::Internal>> mDot12P
         (
-            cavitation_->mDotP()
+            cavitation_->mDot12P()
         );
 
-        const volScalarField::Internal vDotcP(pCoeff*mDotP[0]);
+        const volScalarField::Internal vDot1P(pCoeff*mDot12P[0]);
-        const volScalarField::Internal vDotvP(pCoeff*mDotP[1]);
+        const volScalarField::Internal vDot2P(pCoeff*mDot12P[1]);
 
         eqn +=
-            (vDotvP - vDotcP)*(cavitation_->pSat() - rho*gh)
+            (vDot2P - vDot1P)*(cavitation_->pSat() - rho*gh)
-          - fvm::Sp(vDotvP - vDotcP, eqn.psi());
+          - fvm::Sp(vDot2P - vDot1P, eqn.psi());
     }
 }
 

applications/solvers/multiphase/interFoam/fvModels/VoFCavitation/VoFCavitation.C

@@ -95,21 +95,21 @@ void Foam::fv::VoFCavitation::addSup
 
     if (fieldName == alphaName_)
     {
-        volScalarField::Internal alphalCoeff
+        const volScalarField::Internal alpha1Coeff
         (
             1.0/mixture_.rho1()
           - mixture_.alpha1()()*(1.0/mixture_.rho1() - 1.0/mixture_.rho2())
         );
 
-        const Pair<tmp<volScalarField::Internal>> mDotAlphal
+        const Pair<tmp<volScalarField::Internal>> mDot12Alpha
         (
-            cavitation_->mDotAlphal()
+            cavitation_->mDot12Alpha()
         );
 
-        const volScalarField::Internal vDotcAlphal(alphalCoeff*mDotAlphal[0]());
+        const volScalarField::Internal vDot1Alpha(alpha1Coeff*mDot12Alpha[0]());
-        const volScalarField::Internal vDotvAlphal(alphalCoeff*mDotAlphal[1]());
+        const volScalarField::Internal vDot2Alpha(alpha1Coeff*mDot12Alpha[1]());
 
-        eqn += fvm::Sp(vDotvAlphal - vDotcAlphal, eqn.psi()) + vDotcAlphal;
+        eqn += fvm::Sp(-vDot2Alpha - vDot1Alpha, eqn.psi()) + vDot1Alpha;
     }
 }
 
@@ -139,17 +139,17 @@ void Foam::fv::VoFCavitation::addSup
             1.0/mixture_.rho1() - 1.0/mixture_.rho2()
         );
 
-        const Pair<tmp<volScalarField::Internal>> mDotP
+        const Pair<tmp<volScalarField::Internal>> mDot12P
         (
-            cavitation_->mDotP()
+            cavitation_->mDot12P()
         );
 
-        const volScalarField::Internal vDotcP(pCoeff*mDotP[0]);
+        const volScalarField::Internal vDot1P(pCoeff*mDot12P[0]);
-        const volScalarField::Internal vDotvP(pCoeff*mDotP[1]);
+        const volScalarField::Internal vDot2P(pCoeff*mDot12P[1]);
 
         eqn +=
-            (vDotvP - vDotcP)*(cavitation_->pSat() - rho*gh)
+            (vDot2P - vDot1P)*(cavitation_->pSat() - rho*gh)
-          - fvm::Sp(vDotvP - vDotcP, eqn.psi());
+          - fvm::Sp(vDot2P - vDot1P, eqn.psi());
     }
 }
 

src/twoPhaseModels/compressibleCavitation/Kunz/Kunz.C

@@ -65,55 +65,53 @@ Foam::compressible::cavitationModels::Kunz::Kunz
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::compressible::cavitationModels::Kunz::mDotAlphal() const
+Foam::compressible::cavitationModels::Kunz::mDotcvAlphal() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
-    const volScalarField::Internal mcCoeff_(Cc_*phases_.rho2()()/tInf_);
+    const volScalarField::Internal mcCoeff_(Cc_*rhov()/tInf_);
     const volScalarField::Internal mvCoeff_
     (
-        Cv_*phases_.rho2()()/(0.5*phases_.rho1()()*sqr(UInf_)*tInf_)
+        Cv_*rhov()/(0.5*rhol()*sqr(UInf_)*tInf_)
     );
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        mcCoeff_*sqr(limitedAlpha1)
+        mcCoeff_*sqr(limitedAlphal)
        *max(p - pSat(), p0_)/max(p - pSat(), 0.01*pSat()),
-
+       -mvCoeff_*min(p - pSat(), p0_)
-        mvCoeff_*min(p - pSat(), p0_)
     );
 }
 
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::compressible::cavitationModels::Kunz::mDotP() const
+Foam::compressible::cavitationModels::Kunz::mDotcvP() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
-    const volScalarField::Internal mcCoeff_(Cc_*phases_.rho2()()/tInf_);
+    const volScalarField::Internal mcCoeff_(Cc_*rhov()/tInf_);
     const volScalarField::Internal mvCoeff_
     (
-        Cv_*phases_.rho2()()/(0.5*phases_.rho1()()*sqr(UInf_)*tInf_)
+        Cv_*rhov()/(0.5*rhol()*sqr(UInf_)*tInf_)
     );
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        mcCoeff_*sqr(limitedAlpha1)*(1.0 - limitedAlpha1)
+        mcCoeff_*sqr(limitedAlphal)*(1.0 - limitedAlphal)
        *pos0(p - pSat())/max(p - pSat(), 0.01*pSat()),
-
+        (-mvCoeff_)*limitedAlphal*neg(p - pSat())
-        (-mvCoeff_)*limitedAlpha1*neg(p - pSat())
     );
 }
 

src/twoPhaseModels/compressibleCavitation/Kunz/Kunz.H

@@ -102,13 +102,13 @@ public:
     // Member Functions
 
         //- Return the mass condensation and vaporisation rates as a
-        //  coefficient to multiply (1 - alphal) for the condensation rate
+        //  coefficient to multiply alphav for the condensation rate and a
-        //  and a coefficient to multiply  alphal for the vaporisation rate
+        //  coefficient to multiply alphal for the vaporisation rate
-        virtual Pair<tmp<volScalarField::Internal>> mDotAlphal() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvAlphal() const;
 
         //- Return the mass condensation and vaporisation rates as coefficients
         //  to multiply (p - pSat)
-        virtual Pair<tmp<volScalarField::Internal>> mDotP() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const;
 
         //- Correct the Kunz phaseChange model
         virtual void correct();

src/twoPhaseModels/compressibleCavitation/Merkle/Merkle.C

@@ -67,45 +67,44 @@ Foam::compressible::cavitationModels::Merkle::Merkle
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::compressible::cavitationModels::Merkle::mDotAlphal() const
+Foam::compressible::cavitationModels::Merkle::mDotcvAlphal() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
     const volScalarField::Internal mvCoeff_
     (
-        Cv_*phases_.rho1()()/(0.5*sqr(UInf_)*tInf_*phases_.rho2()())
+        Cv_*rhol()/(0.5*sqr(UInf_)*tInf_*rhov())
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
         mcCoeff_*max(p - pSat(), p0_),
-        mvCoeff_*min(p - pSat(), p0_)
+       -mvCoeff_*min(p - pSat(), p0_)
     );
 }
 
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::compressible::cavitationModels::Merkle::mDotP() const
+Foam::compressible::cavitationModels::Merkle::mDotcvP() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     const volScalarField::Internal mvCoeff_
     (
-        Cv_*phases_.rho1()()/(0.5*sqr(UInf_)*tInf_*phases_.rho2()())
+        Cv_*rhol()/(0.5*sqr(UInf_)*tInf_*rhov())
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        mcCoeff_*(1.0 - limitedAlpha1)*pos0(p - pSat()),
+        mcCoeff_*(1.0 - limitedAlphal)*pos0(p - pSat()),
-
+        (-mvCoeff_)*limitedAlphal*neg(p - pSat())
-        (-mvCoeff_)*limitedAlpha1*neg(p - pSat())
     );
 }
 

src/twoPhaseModels/compressibleCavitation/Merkle/Merkle.H

@@ -98,13 +98,13 @@ public:
     // Member Functions
 
         //- Return the mass condensation and vaporisation rates as a
-        //  coefficient to multiply (1 - alphal) for the condensation rate
+        //  coefficient to multiply alphav for the condensation rate and a
-        //  and a coefficient to multiply  alphal for the vaporisation rate
+        //  coefficient to multiply alphal for the vaporisation rate
-        virtual Pair<tmp<volScalarField::Internal>> mDotAlphal() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvAlphal() const;
 
         //- Return the mass condensation and vaporisation rates as coefficients
         //  to multiply (p - pSat)
-        virtual Pair<tmp<volScalarField::Internal>> mDotP() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const;
 
         //- Correct the Merkle phaseChange model
         virtual void correct();

src/twoPhaseModels/compressibleCavitation/SchnerrSauer/SchnerrSauer.C

@@ -73,13 +73,13 @@ Foam::compressible::cavitationModels::SchnerrSauer::SchnerrSauer
 Foam::tmp<Foam::volScalarField::Internal>
 Foam::compressible::cavitationModels::SchnerrSauer::rRb
 (
-    const volScalarField::Internal& limitedAlpha1
+    const volScalarField::Internal& limitedAlphal
 ) const
 {
     return pow
     (
         ((4*constant::mathematical::pi*n_)/3)
-       *limitedAlpha1/(1.0 + alphaNuc() - limitedAlpha1),
+       *limitedAlphal/(1.0 + alphaNuc() - limitedAlphal),
         1.0/3.0
     );
 }
@@ -99,65 +99,63 @@ Foam::compressible::cavitationModels::SchnerrSauer::pCoeff
     const volScalarField::Internal& p
 ) const
 {
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     const volScalarField::Internal rho
     (
-        limitedAlpha1*phases_.rho1()()
+        limitedAlphal*rhol()
-      + (scalar(1) - limitedAlpha1)*phases_.rho2()()
+      + (scalar(1) - limitedAlphal)*rhov()
     );
 
     return
-        (3*phases_.rho1()()*phases_.rho2()())*sqrt(2/(3*phases_.rho1()()))
+        (3*rhol()*rhov())*sqrt(2/(3*rhol()))
-       *rRb(limitedAlpha1)/(rho*sqrt(mag(p - pSat()) + 0.01*pSat()));
+       *rRb(limitedAlphal)/(rho*sqrt(mag(p - pSat()) + 0.01*pSat()));
 }
 
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::compressible::cavitationModels::SchnerrSauer::mDotAlphal() const
+Foam::compressible::cavitationModels::SchnerrSauer::mDotcvAlphal() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
     const volScalarField::Internal pCoeff(this->pCoeff(p));
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        Cc_*limitedAlpha1*pCoeff*max(p - pSat(), p0_),
+        Cc_*limitedAlphal*pCoeff*max(p - pSat(), p0_),
-
+       -Cv_*(1.0 + alphaNuc() - limitedAlphal)*pCoeff*min(p - pSat(), p0_)
-        Cv_*(1.0 + alphaNuc() - limitedAlpha1)*pCoeff*min(p - pSat(), p0_)
     );
 }
 
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::compressible::cavitationModels::SchnerrSauer::mDotP() const
+Foam::compressible::cavitationModels::SchnerrSauer::mDotcvP() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
     const volScalarField::Internal pCoeff(this->pCoeff(p));
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
-    const volScalarField::Internal apCoeff(limitedAlpha1*pCoeff);
+    const volScalarField::Internal apCoeff(limitedAlphal*pCoeff);
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        Cc_*(1.0 - limitedAlpha1)*pos0(p - pSat())*apCoeff,
+        Cc_*(1.0 - limitedAlphal)*pos0(p - pSat())*apCoeff,
-
+        (-Cv_)*(1.0 + alphaNuc() - limitedAlphal)*neg(p - pSat())*apCoeff
-        (-Cv_)*(1.0 + alphaNuc() - limitedAlpha1)*neg(p - pSat())*apCoeff
     );
 }
 

src/twoPhaseModels/compressibleCavitation/SchnerrSauer/SchnerrSauer.H

@@ -82,9 +82,9 @@ class SchnerrSauer
         dimensionedScalar alphaNuc() const;
 
         //- Reciprocal bubble radius
-        tmp<volScalarField::Internal>rRb
+        tmp<volScalarField::Internal> rRb
         (
-            const volScalarField::Internal& limitedAlpha1
+            const volScalarField::Internal& limitedAlphal
         ) const;
 
         //- Part of the condensation and vapourisation rates
@@ -118,13 +118,13 @@ public:
     // Member Functions
 
         //- Return the mass condensation and vaporisation rates as a
-        //  coefficient to multiply (1 - alphal) for the condensation rate
+        //  coefficient to multiply alphav for the condensation rate and a
-        //  and a coefficient to multiply  alphal for the vaporisation rate
+        //  coefficient to multiply alphal for the vaporisation rate
-        virtual Pair<tmp<volScalarField::Internal>> mDotAlphal() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvAlphal() const;
 
         //- Return the mass condensation and vaporisation rates as coefficients
         //  to multiply (p - pSat)
-        virtual Pair<tmp<volScalarField::Internal>> mDotP() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const;
 
         //- Correct the SchnerrSauer phaseChange model
         virtual void correct();

src/twoPhaseModels/compressibleCavitation/cavitationModel/cavitationModel.C

@@ -47,6 +47,7 @@ Foam::compressible::cavitationModel::cavitationModel
 )
 :
     phases_(phases),
+    liquidIndex_(phases.index(dict.lookup<word>("liquid"))),
     pSat_("pSat", dimPressure, dict.lookup("pSat"))
 {}
 

src/twoPhaseModels/compressibleCavitation/cavitationModel/cavitationModel.H

@@ -56,13 +56,43 @@ protected:
 
     // Protected data
 
-        //- Mixture properties
+        //- Phases
         const compressibleTwoPhases& phases_;
 
+        //- Index of the liquid
+        const bool liquidIndex_;
+
         //- Saturation vapour pressure
         dimensionedScalar pSat_;
 
 
+    // Protected Member Functions
+
+        //- Return the liquid density
+        inline const volScalarField::Internal& alphal() const
+        {
+            return phases_.alpha(liquidIndex_);
+        }
+
+        //- Return the vapour density
+        inline const volScalarField::Internal& alphav() const
+        {
+            return phases_.alpha(!liquidIndex_);
+        }
+
+        //- Return the liquid density
+        inline const volScalarField::Internal& rhol() const
+        {
+            return phases_.rho(liquidIndex_);
+        }
+
+        //- Return the vapour density
+        inline const volScalarField::Internal& rhov() const
+        {
+            return phases_.rho(!liquidIndex_);
+        }
+
+
 public:
 
     //- Runtime type information
@@ -109,20 +139,34 @@ public:
 
     // Member Functions
 
-        //- Return const-access to the saturation vapour pressure
+        //- Return the saturation vapour pressure
-        const dimensionedScalar& pSat() const
+        inline const dimensionedScalar& pSat() const
         {
             return pSat_;
         }
 
         //- Return the mass condensation and vaporisation rates as a
-        //  coefficient to multiply (1 - alphal) for the condensation rate
+        //  coefficient to multiply alphav for the condensation rate and a
-        //  and a coefficient to multiply  alphal for the vaporisation rate
+        //  coefficient to multiply alphal for the vaporisation rate
-        virtual Pair<tmp<volScalarField::Internal>> mDotAlphal() const = 0;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvAlphal() const = 0;
 
         //- Return the mass condensation and vaporisation rates as coefficients
         //  to multiply (p - pSat)
-        virtual Pair<tmp<volScalarField::Internal>> mDotP() const = 0;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const = 0;
+
+        //- Return the mass transfer rates of the two phases as coefficients to
+        //  multiply the volume fraction of the other phase
+        inline Pair<tmp<volScalarField::Internal>> mDot12Alpha() const
+        {
+            return liquidIndex_ ? reverse(mDotcvAlphal()) : mDotcvAlphal();
+        }
+
+        //- Return the mass transfer rates of the two phases as coefficients to
+        //  multiply (p - pSat)
+        inline Pair<tmp<volScalarField::Internal>> mDot12P() const
+        {
+            return liquidIndex_ ? reverse(mDotcvP()) : mDotcvP();
+        }
 
         //- Correct the cavitation model
         virtual void correct() = 0;

src/twoPhaseModels/compressibleTwoPhaseMixture/compressibleTwoPhases.H

@@ -77,11 +77,23 @@ public:
         //- Return the thermo for phase 2
         virtual const rhoThermo& thermo2() const = 0;
 
+        //- Return the density of a given phase
+        inline const rhoThermo& thermo(const bool index) const
+        {
+            return index ? thermo2() : thermo1();
+        }
+
         //- Return the density of phase 1
         virtual const volScalarField& rho1() const = 0;
 
         //- Return the density of phase 2
         virtual const volScalarField& rho2() const = 0;
+
+        //- Return the density of a given phase
+        inline const volScalarField& rho(const bool index) const
+        {
+            return index ? rho2() : rho1();
+        }
 };
 
 

src/twoPhaseModels/incompressibleCavitation/Kunz/Kunz.C

@@ -54,8 +54,8 @@ Foam::cavitationModels::Kunz::Kunz
 
     p0_("0", pSat().dimensions(), 0.0),
 
-    mcCoeff_(Cc_*phases_.rho2()/tInf_),
+    mcCoeff_(Cc_*rhov()/tInf_),
-    mvCoeff_(Cv_*phases_.rho2()/(0.5*phases_.rho1()*sqr(UInf_)*tInf_))
+    mvCoeff_(Cv_*rhov()/(0.5*rhol()*sqr(UInf_)*tInf_))
 {
     correct();
 }
@@ -64,43 +64,41 @@ Foam::cavitationModels::Kunz::Kunz
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::cavitationModels::Kunz::mDotAlphal() const
+Foam::cavitationModels::Kunz::mDotcvAlpha() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        mcCoeff_*sqr(limitedAlpha1)
+        mcCoeff_*sqr(limitedAlphal)
        *max(p - pSat(), p0_)/max(p - pSat(), 0.01*pSat()),
-
+       -mvCoeff_*min(p - pSat(), p0_)
-        mvCoeff_*min(p - pSat(), p0_)
     );
 }
 
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::cavitationModels::Kunz::mDotP() const
+Foam::cavitationModels::Kunz::mDotcvP() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        mcCoeff_*sqr(limitedAlpha1)*(1.0 - limitedAlpha1)
+        mcCoeff_*sqr(limitedAlphal)*(1.0 - limitedAlphal)
        *pos0(p - pSat())/max(p - pSat(), 0.01*pSat()),
-
+        (-mvCoeff_)*limitedAlphal*neg(p - pSat())
-        (-mvCoeff_)*limitedAlpha1*neg(p - pSat())
     );
 }
 
@@ -118,8 +116,8 @@ bool Foam::cavitationModels::Kunz::read(const dictionary& dict)
         dict.lookup("Cc") >> Cc_;
         dict.lookup("Cv") >> Cv_;
 
-        mcCoeff_ = Cc_*phases_.rho2()/tInf_;
+        mcCoeff_ = Cc_*rhov()/tInf_;
-        mvCoeff_ = Cv_*phases_.rho2()/(0.5*phases_.rho1()*sqr(UInf_)*tInf_);
+        mvCoeff_ = Cv_*rhov()/(0.5*rhol()*sqr(UInf_)*tInf_);
 
         return true;
     }

src/twoPhaseModels/incompressibleCavitation/Kunz/Kunz.H

@@ -103,13 +103,13 @@ public:
     // Member Functions
 
         //- Return the mass condensation and vaporisation rates as a
-        //  coefficient to multiply (1 - alphal) for the condensation rate
+        //  coefficient to multiply alphav for the condensation rate and a
-        //  and a coefficient to multiply  alphal for the vaporisation rate
+        //  coefficient to multiply alphal for the vaporisation rate
-        virtual Pair<tmp<volScalarField::Internal>> mDotAlphal() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvAlpha() const;
 
         //- Return the mass condensation and vaporisation rates as coefficients
         //  to multiply (p - pSat)
-        virtual Pair<tmp<volScalarField::Internal>> mDotP() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const;
 
         //- Correct the Kunz phaseChange model
         virtual void correct();

src/twoPhaseModels/incompressibleCavitation/Merkle/Merkle.C

@@ -55,7 +55,7 @@ Foam::cavitationModels::Merkle::Merkle
     p0_("0", pSat().dimensions(), 0.0),
 
     mcCoeff_(Cc_/(0.5*sqr(UInf_)*tInf_)),
-    mvCoeff_(Cv_*phases_.rho1()/(0.5*sqr(UInf_)*tInf_*phases_.rho2()))
+    mvCoeff_(Cv_*rhol()/(0.5*sqr(UInf_)*tInf_*rhov()))
 {
     correct();
 }
@@ -64,7 +64,7 @@ Foam::cavitationModels::Merkle::Merkle
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::cavitationModels::Merkle::mDotAlphal() const
+Foam::cavitationModels::Merkle::mDotcvAlpha() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
@@ -72,26 +72,26 @@ Foam::cavitationModels::Merkle::mDotAlphal() const
     return Pair<tmp<volScalarField::Internal>>
     (
         mcCoeff_*max(p - pSat(), p0_),
-        mvCoeff_*min(p - pSat(), p0_)
+       -mvCoeff_*min(p - pSat(), p0_)
     );
 }
 
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::cavitationModels::Merkle::mDotP() const
+Foam::cavitationModels::Merkle::mDotcvP() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        mcCoeff_*(1.0 - limitedAlpha1)*pos0(p - pSat()),
+        mcCoeff_*(1.0 - limitedAlphal)*pos0(p - pSat()),
-        (-mvCoeff_)*limitedAlpha1*neg(p - pSat())
+        (-mvCoeff_)*limitedAlphal*neg(p - pSat())
     );
 }
 
@@ -110,7 +110,7 @@ bool Foam::cavitationModels::Merkle::read(const dictionary& dict)
         dict.lookup("Cv") >> Cv_;
 
         mcCoeff_ = Cc_/(0.5*sqr(UInf_)*tInf_);
-        mvCoeff_ = Cv_*phases_.rho1()/(0.5*sqr(UInf_)*tInf_*phases_.rho2());
+        mvCoeff_ = Cv_*rhol()/(0.5*sqr(UInf_)*tInf_*rhov());
 
         return true;
     }

src/twoPhaseModels/incompressibleCavitation/Merkle/Merkle.H

@@ -97,13 +97,13 @@ public:
     // Member Functions
 
         //- Return the mass condensation and vaporisation rates as a
-        //  coefficient to multiply (1 - alphal) for the condensation rate
+        //  coefficient to multiply alphav for the condensation rate and a
-        //  and a coefficient to multiply  alphal for the vaporisation rate
+        //  coefficient to multiply alphal for the vaporisation rate
-        virtual Pair<tmp<volScalarField::Internal>> mDotAlphal() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvAlpha() const;
 
         //- Return the mass condensation and vaporisation rates as coefficients
         //  to multiply (p - pSat)
-        virtual Pair<tmp<volScalarField::Internal>> mDotP() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const;
 
         //- Correct the Merkle phaseChange model
         virtual void correct();

src/twoPhaseModels/incompressibleCavitation/SchnerrSauer/SchnerrSauer.C

@@ -70,13 +70,13 @@ Foam::cavitationModels::SchnerrSauer::SchnerrSauer
 Foam::tmp<Foam::volScalarField::Internal>
 Foam::cavitationModels::SchnerrSauer::rRb
 (
-    const volScalarField::Internal& limitedAlpha1
+    const volScalarField::Internal& limitedAlphal
 ) const
 {
     return pow
     (
         ((4*constant::mathematical::pi*n_)/3)
-       *limitedAlpha1/(1.0 + alphaNuc() - limitedAlpha1),
+       *limitedAlphal/(1.0 + alphaNuc() - limitedAlphal),
         1.0/3.0
     );
 }
@@ -96,65 +96,63 @@ Foam::cavitationModels::SchnerrSauer::pCoeff
     const volScalarField::Internal& p
 ) const
 {
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     const volScalarField::Internal rho
     (
-        limitedAlpha1*phases_.rho1()
+        limitedAlphal*rhol()
-      + (scalar(1) - limitedAlpha1)*phases_.rho2()
+      + (scalar(1) - limitedAlphal)*rhov()
     );
 
     return
-        (3*phases_.rho1()*phases_.rho2())*sqrt(2/(3*phases_.rho1()))
+        (3*rhol()*rhov())*sqrt(2/(3*rhol()))
-       *rRb(limitedAlpha1)/(rho*sqrt(mag(p - pSat()) + 0.01*pSat()));
+       *rRb(limitedAlphal)/(rho*sqrt(mag(p - pSat()) + 0.01*pSat()));
 }
 
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::cavitationModels::SchnerrSauer::mDotAlphal() const
+Foam::cavitationModels::SchnerrSauer::mDotcvAlpha() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
     const volScalarField::Internal pCoeff(this->pCoeff(p));
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        Cc_*limitedAlpha1*pCoeff*max(p - pSat(), p0_),
+        Cc_*limitedAlphal*pCoeff*max(p - pSat(), p0_),
-
+       -Cv_*(1.0 + alphaNuc() - limitedAlphal)*pCoeff*min(p - pSat(), p0_)
-        Cv_*(1.0 + alphaNuc() - limitedAlpha1)*pCoeff*min(p - pSat(), p0_)
     );
 }
 
 
 Foam::Pair<Foam::tmp<Foam::volScalarField::Internal>>
-Foam::cavitationModels::SchnerrSauer::mDotP() const
+Foam::cavitationModels::SchnerrSauer::mDotcvP() const
 {
     const volScalarField::Internal& p =
         phases_.mesh().lookupObject<volScalarField>("p");
 
     const volScalarField::Internal pCoeff(this->pCoeff(p));
 
-    const volScalarField::Internal limitedAlpha1
+    const volScalarField::Internal limitedAlphal
     (
-        min(max(phases_.alpha1()(), scalar(0)), scalar(1))
+        min(max(alphal(), scalar(0)), scalar(1))
     );
 
-    const volScalarField::Internal apCoeff(limitedAlpha1*pCoeff);
+    const volScalarField::Internal apCoeff(limitedAlphal*pCoeff);
 
     return Pair<tmp<volScalarField::Internal>>
     (
-        Cc_*(1.0 - limitedAlpha1)*pos0(p - pSat())*apCoeff,
+        Cc_*(1.0 - limitedAlphal)*pos0(p - pSat())*apCoeff,
-
+        (-Cv_)*(1.0 + alphaNuc() - limitedAlphal)*neg(p - pSat())*apCoeff
-        (-Cv_)*(1.0 + alphaNuc() - limitedAlpha1)*neg(p - pSat())*apCoeff
     );
 }
 

src/twoPhaseModels/incompressibleCavitation/SchnerrSauer/SchnerrSauer.H

@@ -80,9 +80,9 @@ class SchnerrSauer
         dimensionedScalar alphaNuc() const;
 
         //- Reciprocal bubble radius
-        tmp<volScalarField::Internal>rRb
+        tmp<volScalarField::Internal> rRb
         (
-            const volScalarField::Internal& limitedAlpha1
+            const volScalarField::Internal& limitedAlphal
         ) const;
 
         //- Part of the condensation and vapourisation rates
@@ -116,13 +116,13 @@ public:
     // Member Functions
 
         //- Return the mass condensation and vaporisation rates as a
-        //  coefficient to multiply (1 - alphal) for the condensation rate
+        //  coefficient to multiply alphav for the condensation rate and a
-        //  and a coefficient to multiply  alphal for the vaporisation rate
+        //  coefficient to multiply alphal for the vaporisation rate
-        virtual Pair<tmp<volScalarField::Internal>> mDotAlphal() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvAlpha() const;
 
         //- Return the mass condensation and vaporisation rates as coefficients
         //  to multiply (p - pSat)
-        virtual Pair<tmp<volScalarField::Internal>> mDotP() const;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const;
 
         //- Correct the SchnerrSauer phaseChange model
         virtual void correct();

src/twoPhaseModels/incompressibleCavitation/cavitationModel/cavitationModel.C

@@ -44,6 +44,7 @@ Foam::cavitationModel::cavitationModel
 )
 :
     phases_(phases),
+    liquidIndex_(phases.index(dict.lookup<word>("liquid"))),
     pSat_("pSat", dimPressure, dict.lookup("pSat"))
 {}
 

src/twoPhaseModels/incompressibleCavitation/cavitationModel/cavitationModel.H

@@ -57,10 +57,40 @@ protected:
         //- Mixture properties
         const incompressibleTwoPhases& phases_;
 
+        //- Index of the liquid
+        const bool liquidIndex_;
+
         //- Saturation vapour pressure
         dimensionedScalar pSat_;
 
 
+    // Protected Member Functions
+
+        //- Return the liquid density
+        inline const volScalarField::Internal& alphal() const
+        {
+            return phases_.alpha(liquidIndex_);
+        }
+
+        //- Return the vapour density
+        inline const volScalarField::Internal& alphav() const
+        {
+            return phases_.alpha(!liquidIndex_);
+        }
+
+        //- Return the liquid density
+        inline const dimensionedScalar& rhol() const
+        {
+            return phases_.rho(liquidIndex_);
+        }
+
+        //- Return the vapour density
+        inline const dimensionedScalar& rhov() const
+        {
+            return phases_.rho(!liquidIndex_);
+        }
+
+
 public:
 
     //- Runtime type information
@@ -107,20 +137,34 @@ public:
 
     // Member Functions
 
-        //- Return const-access to the saturation vapour pressure
+        //- Return the saturation vapour pressure
-        const dimensionedScalar& pSat() const
+        inline const dimensionedScalar& pSat() const
         {
             return pSat_;
         }
 
         //- Return the mass condensation and vaporisation rates as a
-        //  coefficient to multiply (1 - alphal) for the condensation rate
+        //  coefficient to multiply alphav for the condensation rate and a
-        //  and a coefficient to multiply  alphal for the vaporisation rate
+        //  coefficient to multiply alphal for the vaporisation rate
-        virtual Pair<tmp<volScalarField::Internal>> mDotAlphal() const = 0;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvAlpha() const = 0;
 
         //- Return the mass condensation and vaporisation rates as coefficients
         //  to multiply (p - pSat)
-        virtual Pair<tmp<volScalarField::Internal>> mDotP() const = 0;
+        virtual Pair<tmp<volScalarField::Internal>> mDotcvP() const = 0;
+
+        //- Return the mass transfer rates of the two phases as coefficients to
+        //  multiply the volume fraction of the other phase
+        inline Pair<tmp<volScalarField::Internal>> mDot12Alpha() const
+        {
+            return liquidIndex_ ? reverse(mDotcvAlpha()) : mDotcvAlpha();
+        }
+
+        //- Return the mass transfer rates of the two phases as coefficients to
+        //  multiply (p - pSat)
+        inline Pair<tmp<volScalarField::Internal>> mDot12P() const
+        {
+            return liquidIndex_ ? reverse(mDotcvP()) : mDotcvP();
+        }
 
         //- Correct the cavitation model
         virtual void correct() = 0;

src/twoPhaseModels/incompressibleTwoPhaseMixture/incompressibleTwoPhases.H

@@ -52,6 +52,7 @@ class incompressibleTwoPhases
 {
 public:
 
+    //- Runtime type information
     TypeName("incompressibleTwoPhases");
 
 
@@ -69,11 +70,17 @@ public:
 
     // Member Functions
 
-        //- Return const-access to phase1 density
+        //- Return the density of phase1
         virtual const dimensionedScalar& rho1() const = 0;
 
-        //- Return const-access to phase2 density
+        //- Return the density of phase2
         virtual const dimensionedScalar& rho2() const = 0;
+
+        //- Return the density of a given phase
+        inline const dimensionedScalar& rho(const bool index) const
+        {
+            return index ? rho2() : rho1();
+        }
 };
 
 

src/twoPhaseModels/twoPhaseMixture/twoPhases.H

@@ -50,6 +50,7 @@ class twoPhases
 {
 public:
 
+    //- Runtime type information
     TypeName("twoPhases");
 
 
@@ -73,11 +74,37 @@ public:
             return alpha1().mesh();
         }
 
-        //- Return the phase-fraction of phase 1
+        //- Return the index associated with a given phase name
+        inline bool index(const word& name) const
+        {
+            if (alpha1().group() == name)
+            {
+                return false;
+            }
+
+            if (alpha2().group() == name)
+            {
+                return true;
+            }
+
+            FatalErrorInFunction
+                << "Index of phase " << name << " requested from phases"
+                << alpha1().group() << " and " << alpha2().group()
+                << exit(FatalError);
+            return false;
+        }
+
+        //- Return the volume-fraction of phase 1
         virtual const volScalarField& alpha1() const = 0;
 
-        //- Return the phase-fraction of phase 2
+        //- Return the volume-fraction of phase 2
         virtual const volScalarField& alpha2() const = 0;
+
+        //- Return the volume-fraction of a given phase
+        const volScalarField& alpha(const bool index) const
+        {
+            return index ? alpha2() : alpha1();
+        }
 };
 
 

tutorials/modules/compressibleVoF/ballValve/0/alpha.vapour

@@ -10,13 +10,13 @@ FoamFile
     format      ascii;
     class       volScalarField;
     location    "0";
-    object      alpha.water;
+    object      alpha.vapour;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 dimensions      [0 0 0 0 0 0 0];
 
-internalField   uniform 1;
+internalField   uniform 0;
 
 boundaryField
 {

tutorials/modules/compressibleVoF/ballValve/constant/fvModels

@@ -22,6 +22,8 @@ VoFCavitation
 
     model   SchnerrSauer;
 
+    liquid  water;
+
     pSat    2300;
 
     n       1.6e+13;

tutorials/modules/compressibleVoF/ballValve/constant/phaseProperties

@@ -14,7 +14,7 @@ FoamFile
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-phases          (water vapour);
+phases          (vapour water);
 
 pMin            100;
 

tutorials/modules/compressibleVoF/ballValve/system/fvSolution

@@ -16,7 +16,7 @@ FoamFile
 
 solvers
 {
-    "alpha.water.*"
+    "alpha.vapour.*"
     {
         nAlphaCorr      2;
         nAlphaSubCycles 1;

tutorials/multiphase/interFoam/RAS/propeller/constant/fvModels

@@ -26,6 +26,8 @@ VoFCavitation
     {
         pSat    2300;   // Saturation pressure
 
+        liquid  water;
+
         UInf    20.0;
         tInf    0.005; // L = 0.1 m
         Cc      1000;
@@ -36,6 +38,8 @@ VoFCavitation
     {
         pSat    2300;   // Saturation pressure
 
+        liquid  water;
+
         UInf    20.0;
         tInf    0.005;  // L = 0.1 m
         Cc      80;
@@ -46,6 +50,8 @@ VoFCavitation
     {
         pSat    2300;   // Saturation pressure
 
+        liquid  water;
+
         n       1.6e+13;
         dNuc    2.0e-06;
         Cc      1;

tutorials/multiphase/interFoam/laminar/cavitatingBullet/constant/fvModels

@@ -26,6 +26,8 @@ VoFCavitation
     {
         pSat    2300;   // Saturation pressure
 
+        liquid  water;
+
         UInf    20.0;
         tInf    0.005; // L = 0.1 m
         Cc      1000;
@@ -36,6 +38,8 @@ VoFCavitation
     {
         pSat    2300;   // Saturation pressure
 
+        liquid  water;
+
         UInf    20.0;
         tInf    0.005;  // L = 0.1 m
         Cc      80;
@@ -46,6 +50,8 @@ VoFCavitation
     {
         pSat    2300;   // Saturation pressure
 
+        liquid  water;
+
         n       1.6e+13;
         dNuc    2.0e-06;
         Cc      1;