reactingEulerFoam: Multiphase thermal phase change and support for multiple mass transfer mechanisms

- Thermal phase change and wall boiling functionality has been generalized to
  support two- and multi- phase simulations.
- Thermal phase change now also allows purePhaseModel, which simplifies case setup.
- The phaseSystem templates have been restructured in preparation of multiple
  simultaneous mass transfer mechanisms. For example, combination of thermal phase
  and inhomogeneous population balance models.

Patch contributed by VTT Technical Research Centre of Finland Ltd and Institute
of Fluid Dynamics, Helmholtz-Zentrum Dresden - Rossendorf (HZDR).

applications/solvers/multiphase/reactingEulerFoam/Allwclean

@@ -4,6 +4,7 @@ cd ${0%/*} || exit 1    # Run from this directory
 wclean libso phaseSystems
 wclean libso interfacialModels
 wclean libso interfacialCompositionModels
+wclean libso derivedFvPatchFields
 reactingTwoPhaseEulerFoam/Allwclean
 reactingMultiphaseEulerFoam/Allwclean
 

applications/solvers/multiphase/reactingEulerFoam/Allwmake

@@ -9,6 +9,9 @@ wmakeLnInclude interfacialCompositionModels
 wmake $targetType phaseSystems
 wmake $targetType interfacialModels
 wmake $targetType interfacialCompositionModels
+wmake $targetType derivedFvPatchFields
+wmakeLnInclude reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels
+wmakeLnInclude reactingMultiphaseEulerFoam/multiphaseCompressibleTurbulenceModels
 reactingTwoPhaseEulerFoam/Allwmake $targetType $*
 reactingMultiphaseEulerFoam/Allwmake $targetType $*
 

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/Make/files

@@ -0,0 +1,29 @@
+wallBoilingSubModels/partitioningModels/partitioningModel/partitioningModel.C
+wallBoilingSubModels/partitioningModels/partitioningModel/newPartitioningModel.C
+wallBoilingSubModels/partitioningModels/phaseFraction/phaseFraction.C
+wallBoilingSubModels/partitioningModels/Lavieville/Lavieville.C
+wallBoilingSubModels/partitioningModels/cosine/cosine.C
+wallBoilingSubModels/partitioningModels/linear/linear.C
+
+wallBoilingSubModels/nucleationSiteModels/nucleationSiteModel/nucleationSiteModel.C
+wallBoilingSubModels/nucleationSiteModels/nucleationSiteModel/newNucleationSiteModel.C
+wallBoilingSubModels/nucleationSiteModels/LemmertChawla/LemmertChawla.C
+
+wallBoilingSubModels/departureDiameterModels/departureDiameterModel/departureDiameterModel.C
+wallBoilingSubModels/departureDiameterModels/departureDiameterModel/newDepartureDiameterModel.C
+wallBoilingSubModels/departureDiameterModels/TolubinskiKostanchuk/TolubinskiKostanchuk.C
+wallBoilingSubModels/departureDiameterModels/KocamustafaogullariIshii/KocamustafaogullariIshii.C
+
+wallBoilingSubModels/departureFrequencyModels/departureFrequencyModel/departureFrequencyModel.C
+wallBoilingSubModels/departureFrequencyModels/departureFrequencyModel/newDepartureFrequencyModel.C
+wallBoilingSubModels/departureFrequencyModels/Cole/Cole.C
+
+alphatPhaseChangeWallFunction/alphatPhaseChangeWallFunctionFvPatchScalarField.C
+alphatPhaseChangeJayatillekeWallFunction/alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField.C
+alphatFixedDmdtWallBoilingWallFunction/alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField.C
+alphatWallBoilingWallFunction/alphatWallBoilingWallFunctionFvPatchScalarField.C
+copiedFixedValue/copiedFixedValueFvPatchScalarField.C
+fixedMultiPhaseHeatFlux/fixedMultiPhaseHeatFluxFvPatchScalarField.C
+
+
+LIB = $(FOAM_LIBBIN)/libreactingEulerianFvPatchFields

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/Make/options

@@ -0,0 +1,20 @@
+EXE_INC = \
+    -I../phaseSystems/lnInclude \
+    -I../interfacialModels/lnInclude\
+    -I../interfacialCompositionModels/lnInclude \
+    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
+    -I$(LIB_SRC)/transportModels/incompressible/transportModel \
+    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
+    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
+    -I$(LIB_SRC)/TurbulenceModels/phaseCompressible/lnInclude \
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/meshTools/lnInclude
+
+LIB_LIBS = \
+    -lreactingPhaseSystem \
+    -lreactingEulerianInterfacialModels \
+    -lreactingEulerianInterfacialCompositionModels \
+    -lfiniteVolume \
+    -lfvOptions \
+    -lmeshTools

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/alphatFixedDmdtWallBoilingWallFunction/alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -44,6 +44,7 @@ alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField
 )
 :
     alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField(p, iF),
+    vaporPhaseName_("vapor"),
     relax_(1.0),
     fixedDmdt_(0.0),
     L_(0.0)
@@ -61,6 +62,7 @@ alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField
 )
 :
     alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField(p, iF, dict),
+    vaporPhaseName_(dict.lookup("vaporPhase")),
     relax_(dict.lookupOrDefault<scalar>("relax", 1.0)),
     fixedDmdt_(dict.lookupOrDefault<scalar>("fixedDmdt", 0.0)),
     L_(dict.lookupOrDefault<scalar>("L", 0.0))
@@ -116,6 +118,53 @@ alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+bool alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField::
+activePhasePair(const phasePairKey& phasePair) const
+{
+    if (phasePair == phasePairKey(vaporPhaseName_,internalField().group()))
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+const scalarField& alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField::
+dmdt(const phasePairKey& phasePair) const
+{
+    if (activePhasePair(phasePair))
+    {
+        return dmdt_;
+    }
+    else
+    {
+        FatalErrorInFunction
+            << " dmdt requested for invalid phasePair!"
+            << abort(FatalError);
+
+        return mDotL_;
+    }
+}
+
+const scalarField& alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField::
+mDotL(const phasePairKey& phasePair) const
+{
+    if (activePhasePair(phasePair))
+    {
+        return mDotL_;
+    }
+    else
+    {
+        FatalErrorInFunction
+            << " mDotL requested for invalid phasePair!"
+            << abort(FatalError);
+
+        return mDotL_;
+    }
+}
+
 void alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
 {
     if (updated())
@@ -138,6 +187,8 @@ void alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField::write
 ) const
 {
     fvPatchField<scalar>::write(os);
+    os.writeKeyword("vaporPhase") << vaporPhaseName_ << token::END_STATEMENT
+    << nl;
     os.writeKeyword("relax") << relax_ << token::END_STATEMENT << nl;
     os.writeKeyword("fixedDmdt") << fixedDmdt_ << token::END_STATEMENT << nl;
     os.writeKeyword("L") << L_ << token::END_STATEMENT << nl;

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/alphatFixedDmdtWallBoilingWallFunction/alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField.H

@@ -63,6 +63,9 @@ class alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField
 {
     // Private data
 
+        //- name on the phase
+        word vaporPhaseName_;
+
         //- dmdt relaxationFactor
         scalar relax_;
 
@@ -151,6 +154,15 @@ public:
 
     // Member functions
 
+        //- Is there phase change mass transfer for this phasePair
+        virtual bool activePhasePair(const phasePairKey&) const;
+
+        //- Return the rate of phase-change for specific phase pair
+        virtual const scalarField& dmdt(const phasePairKey&) const;
+
+        //- Return the rate of phase-change for specific phase pair
+        virtual const scalarField& mDotL(const phasePairKey&) const;
+
         // Evaluation functions
 
             //- Update the coefficients associated with the patch field

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/alphatPhaseChangeJayatillekeWallFunction/alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -27,9 +27,7 @@ License
 #include "fvPatchFieldMapper.H"
 #include "addToRunTimeSelectionTable.H"
 
-#include "twoPhaseSystem.H"
+#include "phaseSystem.H"
-#include "ThermalPhaseChangePhaseSystem.H"
-#include "MomentumTransferPhaseSystem.H"
 #include "compressibleTurbulenceModel.H"
 #include "ThermalDiffusivity.H"
 #include "PhaseCompressibleTurbulenceModel.H"
@@ -121,33 +119,20 @@ alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField::calcAlphat
     const scalarField& prevAlphat
 ) const
 {
-    // Lookup the fluid model
-    const ThermalPhaseChangePhaseSystem
-    <
-        MomentumTransferPhaseSystem<twoPhaseSystem>
-    >& fluid =
-        refCast
-        <
-            const ThermalPhaseChangePhaseSystem
-            <
-                MomentumTransferPhaseSystem<twoPhaseSystem>
-            >
-        >
-        (
-            db().lookupObject<phaseSystem>("phaseProperties")
-        );
 
-    const phaseModel& liquid
+    // Lookup the fluid model
+    const phaseSystem& fluid =
+        db().lookupObject<phaseSystem>("phaseProperties");
+
+    const phaseModel& phase
     (
-        fluid.phase1().name() == internalField().group()
+        fluid.phases()[internalField().group()]
-      ? fluid.phase1()
-      : fluid.phase2()
     );
 
     const label patchi = patch().index();
 
     // Retrieve turbulence properties from model
-    const phaseCompressibleTurbulenceModel& turbModel = liquid.turbulence();
+    const phaseCompressibleTurbulenceModel& turbModel = phase.turbulence();
 
     const scalar Cmu25 = pow025(Cmu_);
 
@@ -156,7 +141,7 @@ alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField::calcAlphat
     const tmp<scalarField> tmuw = turbModel.mu(patchi);
     const scalarField& muw = tmuw();
 
-    const tmp<scalarField> talphaw = liquid.thermo().alpha(patchi);
+    const tmp<scalarField> talphaw = phase.thermo().alpha(patchi);
     const scalarField& alphaw = talphaw();
 
     const tmp<volScalarField> tk = turbModel.k();
@@ -169,10 +154,10 @@ alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField::calcAlphat
 
     const fvPatchScalarField& rhow = turbModel.rho().boundaryField()[patchi];
     const fvPatchScalarField& hew =
-        liquid.thermo().he().boundaryField()[patchi];
+        phase.thermo().he().boundaryField()[patchi];
 
     const fvPatchScalarField& Tw =
-        liquid.thermo().T().boundaryField()[patchi];
+        phase.thermo().T().boundaryField()[patchi];
 
     scalarField Tp(Tw.patchInternalField());
 

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/alphatPhaseChangeJayatillekeWallFunction/alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/alphatPhaseChangeWallFunction/alphatPhaseChangeWallFunctionFvPatchScalarField.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -120,7 +120,8 @@ alphatPhaseChangeWallFunctionFvPatchScalarField
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-void alphatPhaseChangeWallFunctionFvPatchScalarField::write(Ostream& os) const
+void alphatPhaseChangeWallFunctionFvPatchScalarField::
+write(Ostream& os) const
 {
     fvPatchField<scalar>::write(os);
     dmdt_.writeEntry("dmdt", os);

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/alphatPhaseChangeWallFunction/alphatPhaseChangeWallFunctionFvPatchScalarField.H

@@ -43,6 +43,7 @@ SourceFiles
 #define alphatPhaseChangeWallFunctionFvPatchScalarField_H
 
 #include "fixedValueFvPatchFields.H"
+#include "phasePairKey.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -121,13 +122,43 @@ public:
     // Member functions
 
         //- Return the rate of phase-change
-        const scalarField& dmdt() const
+        virtual const scalarField& dmdt() const
         {
             return dmdt_;
         }
 
-        //- Return the enthelpy source due to phase-change
+        //- Return the enthalpy source due to phase-change
-        const scalarField& mDotL() const
+        virtual const scalarField& mDotL() const
+        {
+            return mDotL_;
+        }
+
+        //- Is there phase change mass transfer for this phasePair
+        virtual bool activePhasePair(const phasePairKey&) const
+        {
+            return false;
+        }
+
+        //- Return the rate of phase-change for specific phase pair
+        virtual const scalarField& dmdt(const phasePairKey&) const
+        {
+            return dmdt_;
+        }
+
+        //- Return the rate of phase-change for specific phase pair
+        virtual const scalarField& mDotL(const phasePairKey&) const
+        {
+            return mDotL_;
+        }
+
+        //- Return the rate of phase-change for specific phase
+        virtual scalarField dmdt(const word&) const
+        {
+            return dmdt_;
+        }
+
+        //- Return the enthalpy source due to phase-change for specific phase
+        virtual scalarField mDotL(const word&) const
         {
             return mDotL_;
         }

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/alphatWallBoilingWallFunction/alphatWallBoilingWallFunctionFvPatchScalarField.C

@@ -27,10 +27,7 @@ License
 #include "fvPatchFieldMapper.H"
 #include "addToRunTimeSelectionTable.H"
 
-#include "twoPhaseSystem.H"
 #include "phaseSystem.H"
-#include "ThermalPhaseChangePhaseSystem.H"
-#include "MomentumTransferPhaseSystem.H"
 #include "compressibleTurbulenceModel.H"
 #include "ThermalDiffusivity.H"
 #include "PhaseCompressibleTurbulenceModel.H"
@@ -82,6 +79,7 @@ alphatWallBoilingWallFunctionFvPatchScalarField
 )
 :
     alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField(p, iF),
+    otherPhaseName_("vapor"),
     phaseType_(liquidPhase),
     relax_(0.5),
     AbyV_(p.size(), 0),
@@ -111,6 +109,7 @@ alphatWallBoilingWallFunctionFvPatchScalarField
 )
 :
     alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField(p, iF, dict),
+    otherPhaseName_(dict.lookup("otherPhase")),
     phaseType_(phaseTypeNames_.read(dict.lookup("phaseType"))),
     relax_(dict.lookupOrDefault<scalar>("relax", 0.5)),
     AbyV_(p.size(), 0),
@@ -122,6 +121,18 @@ alphatWallBoilingWallFunctionFvPatchScalarField
     departureDiamModel_(nullptr),
     departureFreqModel_(nullptr)
 {
+
+    // Check that otherPhaseName != this phase
+    if (internalField().group() == otherPhaseName_)
+    {
+        FatalErrorInFunction
+            << "otherPhase should be the name of the vapor phase that "
+            << "corresponds to the liquid base of vice versa" << nl
+            << "This phase: " << internalField().group() << nl
+            << "otherPhase: " << otherPhaseName_
+            << abort(FatalError);
+    }
+
     switch (phaseType_)
     {
         case vaporPhase:
@@ -206,6 +217,8 @@ alphatWallBoilingWallFunctionFvPatchScalarField
         iF,
         mapper
     ),
+    otherPhaseName_(psf.otherPhaseName_),
+    phaseType_(psf.phaseType_),
     relax_(psf.relax_),
     AbyV_(psf.AbyV_),
     alphatConv_(psf.alphatConv_, mapper),
@@ -225,6 +238,8 @@ alphatWallBoilingWallFunctionFvPatchScalarField
 )
 :
     alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField(psf),
+    otherPhaseName_(psf.otherPhaseName_),
+    phaseType_(psf.phaseType_),
     relax_(psf.relax_),
     AbyV_(psf.AbyV_),
     alphatConv_(psf.alphatConv_),
@@ -245,6 +260,8 @@ alphatWallBoilingWallFunctionFvPatchScalarField
 )
 :
     alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField(psf, iF),
+    otherPhaseName_(psf.otherPhaseName_),
+    phaseType_(psf.phaseType_),
     relax_(psf.relax_),
     AbyV_(psf.AbyV_),
     alphatConv_(psf.alphatConv_),
@@ -259,6 +276,53 @@ alphatWallBoilingWallFunctionFvPatchScalarField
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+bool alphatWallBoilingWallFunctionFvPatchScalarField::
+activePhasePair(const phasePairKey& phasePair) const
+{
+    if (phasePair == phasePairKey(otherPhaseName_,internalField().group()))
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+const scalarField& alphatWallBoilingWallFunctionFvPatchScalarField::
+dmdt(const phasePairKey& phasePair) const
+{
+    if (activePhasePair(phasePair))
+    {
+        return dmdt_;
+    }
+    else
+    {
+        FatalErrorInFunction
+            << " dmdt requested for invalid phasePair!"
+            << abort(FatalError);
+
+        return dmdt_;
+    }
+}
+
+const scalarField& alphatWallBoilingWallFunctionFvPatchScalarField::
+mDotL(const phasePairKey& phasePair) const
+{
+    if (activePhasePair(phasePair))
+    {
+        return mDotL_;
+    }
+    else
+    {
+        FatalErrorInFunction
+            << " mDotL requested for invalid phasePair!"
+            << abort(FatalError);
+
+        return mDotL_;
+    }
+}
+
 void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
 {
     if (updated())
@@ -275,21 +339,15 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
     }
 
     // Lookup the fluid model
-    const ThermalPhaseChangePhaseSystem
+    const phaseSystem& fluid =
-    <
+        refCast<const phaseSystem>
-        MomentumTransferPhaseSystem<twoPhaseSystem>
-    >& fluid =
-        refCast
-        <
-            const ThermalPhaseChangePhaseSystem
-            <
-                MomentumTransferPhaseSystem<twoPhaseSystem>
-            >
-        >
         (
             db().lookupObject<phaseSystem>("phaseProperties")
         );
 
+    const saturationModel& satModel =
+        db().lookupObject<saturationModel>("saturationModel");
+
     const label patchi = patch().index();
 
     switch (phaseType_)
@@ -298,22 +356,17 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
         {
             const phaseModel& vapor
             (
-                fluid.phase1().name() == internalField().group()
+                fluid.phases()[internalField().group()]
-              ? fluid.phase1()
-              : fluid.phase2()
             );
 
-            const phaseModel& liquid(fluid.otherPhase(vapor));
-
-            // Liquid phase fraction at the wall
-            const scalarField liquidw(liquid.boundaryField()[patchi]);
-
             // Vapor Liquid phase fraction at the wall
             const scalarField vaporw(vapor.boundaryField()[patchi]);
 
+            // NOTE! Assumes 1-thisPhase for liquid fraction in
+            // multiphase simulations
             const scalarField fLiquid
             (
-                partitioningModel_->fLiquid(liquidw)
+                partitioningModel_->fLiquid(1-vaporw)
             );
 
             operator==
@@ -350,12 +403,10 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
 
             const phaseModel& liquid
             (
-                fluid.phase1().name() == internalField().group()
+                fluid.phases()[internalField().group()]
-              ? fluid.phase1()
-              : fluid.phase2()
             );
 
-            const phaseModel& vapor(fluid.otherPhase(liquid));
+            const phaseModel& vapor(fluid.phases()[otherPhaseName_]);
 
             // Retrieve turbulence properties from model
             const phaseCompressibleTurbulenceModel& turbModel =
@@ -417,7 +468,8 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
 
             // Saturation temperature
             const tmp<volScalarField> tTsat =
-                fluid.saturation().Tsat(liquid.thermo().p());
+                satModel.Tsat(liquid.thermo().p());
+
             const volScalarField& Tsat = tTsat();
             const fvPatchScalarField& Tsatw(Tsat.boundaryField()[patchi]);
             const scalarField Tsatc(Tsatw.patchInternalField());
@@ -538,16 +590,6 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
 
                 scalar maxErr(max(mag(TsupPrev - TsupNew)));
 
-                if (maxErr < 1e-1)
-                {
-                    if (i > 0)
-                    {
-                        Info<< "Wall boiling wall function iterations: "
-                            << i + 1 << endl;
-                    }
-                    break;
-                }
-
                 if (debug)
                 {
                     const scalarField qc
@@ -586,6 +628,17 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::updateCoeffs()
                     Info<< "  alphatConv: " << gMin(alphatConv_)
                         << " - " << gMax(alphatConv_) << endl;
                 }
+
+                if (maxErr < 1e-1)
+                {
+                    if (i > 0)
+                    {
+                        Info<< "Wall boiling wall function iterations: "
+                            << i + 1 << endl;
+                    }
+                    break;
+                }
+
             }
             break;
         }
@@ -646,6 +699,8 @@ void alphatWallBoilingWallFunctionFvPatchScalarField::write(Ostream& os) const
         }
     }
 
+    os.writeKeyword("otherPhase") << otherPhaseName_ << token::END_STATEMENT
+    << nl;
     dmdt_.writeEntry("dmdt", os);
     dDep_.writeEntry("dDep", os);
     qq_.writeEntry("qQuenching", os);

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/alphatWallBoilingWallFunction/alphatWallBoilingWallFunctionFvPatchScalarField.H

@@ -163,6 +163,9 @@ private:
 
     // Private data
 
+        //- name of the other phase (vapor/liquid phase)
+        word otherPhaseName_;
+
         //- Heat source type names
         static const NamedEnum<phaseType, 2> phaseTypeNames_;
 
@@ -272,6 +275,17 @@ public:
 
     // Member functions
 
+        using alphatPhaseChangeWallFunctionFvPatchScalarField::dmdt;
+
+        //- Is there phase change mass transfer for this phasePair
+        virtual bool activePhasePair(const phasePairKey&) const;
+
+        //- Return the rate of phase-change for specific phase pair
+        virtual const scalarField& dmdt(const phasePairKey&) const;
+
+        //- Return the rate of phase-change for specific phase pair
+        virtual const scalarField& mDotL(const phasePairKey&) const;
+
         //- Return the departure diameter field
         const scalarField& dDeparture() const
         {

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/copiedFixedValue/copiedFixedValueFvPatchScalarField.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/copiedFixedValue/copiedFixedValueFvPatchScalarField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/fixedMultiPhaseHeatFlux/fixedMultiPhaseHeatFluxFvPatchScalarField.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -27,9 +27,7 @@ License
 #include "fvPatchFieldMapper.H"
 #include "addToRunTimeSelectionTable.H"
 
-#include "twoPhaseSystem.H"
+#include "phaseSystem.H"
-#include "ThermalPhaseChangePhaseSystem.H"
-#include "MomentumTransferPhaseSystem.H"
 #include "compressibleTurbulenceModel.H"
 #include "ThermalDiffusivity.H"
 #include "PhaseCompressibleTurbulenceModel.H"
@@ -45,7 +43,8 @@ fixedMultiPhaseHeatFluxFvPatchScalarField
 :
     fixedValueFvPatchScalarField(p, iF),
     q_(p.size(), 0.0),
-    relax_(1.0)
+    relax_(1.0),
+    Tmin_(0.0)
 {}
 
 
@@ -59,7 +58,8 @@ fixedMultiPhaseHeatFluxFvPatchScalarField
 :
     fixedValueFvPatchScalarField(p, iF, dict),
     q_("q", dict, p.size()),
-    relax_(dict.lookupOrDefault<scalar>("relax", 1.0))
+    relax_(dict.lookupOrDefault<scalar>("relax", 1.0)),
+    Tmin_(dict.lookupOrDefault<scalar>("Tmin", 273))
 {}
 
 
@@ -74,7 +74,8 @@ fixedMultiPhaseHeatFluxFvPatchScalarField
 :
     fixedValueFvPatchScalarField(psf, p, iF, mapper),
     q_(psf.q_, mapper),
-    relax_(psf.relax_)
+    relax_(psf.relax_),
+    Tmin_(psf.Tmin_)
 {}
 
 
@@ -86,7 +87,8 @@ fixedMultiPhaseHeatFluxFvPatchScalarField
 :
     fixedValueFvPatchScalarField(psf),
     q_(psf.q_),
-    relax_(psf.relax_)
+    relax_(psf.relax_),
+    Tmin_(psf.Tmin_)
 {}
 
 
@@ -99,7 +101,8 @@ fixedMultiPhaseHeatFluxFvPatchScalarField
 :
     fixedValueFvPatchScalarField(psf, iF),
     q_(psf.q_),
-    relax_(psf.relax_)
+    relax_(psf.relax_),
+    Tmin_(psf.Tmin_)
 {}
 
 
@@ -114,17 +117,7 @@ void Foam::fixedMultiPhaseHeatFluxFvPatchScalarField::updateCoeffs()
     }
 
     // Lookup the fluid model
-    const ThermalPhaseChangePhaseSystem
+    const phaseSystem& fluid =
-    <
-        MomentumTransferPhaseSystem<twoPhaseSystem>
-    >& fluid =
-        refCast
-        <
-            const ThermalPhaseChangePhaseSystem
-            <
-                MomentumTransferPhaseSystem<twoPhaseSystem>
-            >
-        >
         (
             db().lookupObject<phaseSystem>("phaseProperties")
         );
@@ -171,10 +164,11 @@ void Foam::fixedMultiPhaseHeatFluxFvPatchScalarField::updateCoeffs()
     if (debug)
     {
         Info<< patch().name() << " " << ": overall heat flux "
-            << gMin(Q) << " - " << gMax(Q) << endl;
+            << gMin(Q) << " - " << gMax(Q) << " W/m2, power: "
+            << gSum(patch().magSf()*Q) << " W" << endl;
     }
 
-    operator==((1 - relax_)*Tp + relax_*(q_ + A)/(B));
+    operator==((1 - relax_)*Tp + relax_*max(Tmin_,(q_ + A)/(B)));
 
     fixedValueFvPatchScalarField::updateCoeffs();
 }

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/fixedMultiPhaseHeatFlux/fixedMultiPhaseHeatFluxFvPatchScalarField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -70,6 +70,9 @@ class fixedMultiPhaseHeatFluxFvPatchScalarField
         //- Relaxation factor
         scalar relax_;
 
+        //- Minimum temperature limit [K]
+        scalar Tmin_;
+
 
 public:
 

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/IATEsource/wallBoiling.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/departureDiameterModels/departureDiameterModel/departureDiameterModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/departureDiameterModels/departureDiameterModel/newDepartureDiameterModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/Cole/Cole.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/Cole/Cole.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/departureFrequencyModel/departureFrequencyModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/departureFrequencyModel/departureFrequencyModel.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/departureFrequencyModel/newDepartureFrequencyModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/nucleationSiteModels/nucleationSiteModel/newNucleationSiteModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/nucleationSiteModels/nucleationSiteModel/nucleationSiteModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/Lavieville/Lavieville.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/cosine/cosine.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/linear/linear.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/linear/linear.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/partitioningModel/newPartitioningModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/partitioningModel/partitioningModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/partitioningModel/partitioningModel.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/phaseFraction/phaseFraction.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/derivedFvPatchFields/wallBoilingSubModels/partitioningModels/phaseFraction/phaseFraction.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/interfaceCompositionModels/NonRandomTwoLiquid/NonRandomTwoLiquid.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -104,14 +104,16 @@ NonRandomTwoLiquid
     (
         saturationModel::New
         (
-            dict.subDict(species1Name_).subDict("interaction")
+            dict.subDict(species1Name_).subDict("interaction"),
+            pair.phase1().mesh()
         ).ptr()
     );
     saturationModel21_.reset
     (
         saturationModel::New
         (
-            dict.subDict(species2Name_).subDict("interaction")
+            dict.subDict(species2Name_).subDict("interaction"),
+            pair.phase1().mesh()
         ).ptr()
     );
 

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/interfaceCompositionModels/Saturated/Saturated.C

@@ -62,7 +62,8 @@ Foam::interfaceCompositionModels::Saturated<Thermo, OtherThermo>::Saturated
     (
         saturationModel::New
         (
-            dict.subDict("saturationPressure")
+            dict.subDict("saturationPressure"),
+            pair.phase1().mesh()
         )
     )
 {

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/Antoine/Antoine.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -40,9 +40,13 @@ namespace saturationModels
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::saturationModels::Antoine::Antoine(const dictionary& dict)
+Foam::saturationModels::Antoine::Antoine
+(
+    const dictionary& dict,
+    const objectRegistry& db
+)
 :
-    saturationModel(),
+    saturationModel(db),
     A_("A", dimless, dict),
     B_("B", dimTemperature, dict),
     C_("C", dimTemperature, dict)

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/Antoine/Antoine.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -81,7 +81,7 @@ public:
     // Constructors
 
         //- Construct from a dictionary
-        Antoine(const dictionary& dict);
+        Antoine(const dictionary& dict, const objectRegistry& db);
 
 
     //- Destructor

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/AntoineExtended/AntoineExtended.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -47,10 +47,11 @@ namespace saturationModels
 
 Foam::saturationModels::AntoineExtended::AntoineExtended
 (
-    const dictionary& dict
+    const dictionary& dict,
+    const objectRegistry& db
 )
 :
-    Antoine(dict),
+    Antoine(dict, db),
     D_("D", dimless, dict),
     F_("F", dimless, dict),
     E_("E", dimless/pow(dimTemperature, F_), dict)

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/AntoineExtended/AntoineExtended.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -81,7 +81,7 @@ public:
     // Constructors
 
         //- Construct from a dictionary
-        AntoineExtended(const dictionary& dict);
+        AntoineExtended(const dictionary& dict, const objectRegistry& db);
 
 
     //- Destructor

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/ArdenBuck/ArdenBuck.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -57,9 +57,13 @@ Foam::saturationModels::ArdenBuck::xByTC
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::saturationModels::ArdenBuck::ArdenBuck(const dictionary& dict)
+Foam::saturationModels::ArdenBuck::ArdenBuck
+(
+    const dictionary& dict,
+    const objectRegistry& db
+)
 :
-    saturationModel()
+    saturationModel(db)
 {}
 
 

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/ArdenBuck/ArdenBuck.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -66,7 +66,7 @@ public:
     // Constructors
 
         //- Construct from a dictionary
-        ArdenBuck(const dictionary& dict);
+        ArdenBuck(const dictionary& dict, const objectRegistry& db);
 
 
     //- Destructor

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/constantSaturationConditions/constantSaturationConditions.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -46,9 +46,13 @@ namespace saturationModels
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::saturationModels::constantSaturationConditions::
-constantSaturationConditions(const dictionary& dict)
+constantSaturationConditions
+(
+    const dictionary& dict,
+    const objectRegistry& db
+)
 :
-    saturationModel(),
+    saturationModel(db),
     pSat_("pSat", dimPressure, dict),
     Tsat_("Tsat", dimTemperature, dict)
 {}

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/constantSaturationConditions/constantSaturationConditions.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -71,7 +71,11 @@ public:
     // Constructors
 
         //- Construct from a dictionary
-        constantSaturationConditions(const dictionary& dict);
+        constantSaturationConditions
+        (
+            const dictionary& dict,
+            const objectRegistry& db
+        );
 
 
     //- Destructor

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/function1/function1.C

@@ -40,9 +40,13 @@ namespace saturationModels
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::saturationModels::function1::function1(const dictionary& dict)
+Foam::saturationModels::function1::function1
+(
+    const dictionary& dict,
+    const objectRegistry& db
+)
 :
-    saturationModel(),
+    saturationModel(db),
     function_
     (
         Function1<scalar>::New("function", dict)

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/function1/function1.H

@@ -107,7 +107,7 @@ public:
     // Constructors
 
         //- Construct from a dictionary
-        function1(const dictionary& dict);
+        function1(const dictionary& dict, const objectRegistry& db);
 
 
     //- Destructor

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/polynomial/polynomial.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -40,9 +40,13 @@ namespace saturationModels
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::saturationModels::polynomial::polynomial(const dictionary& dict)
+Foam::saturationModels::polynomial::polynomial
+(
+    const dictionary& dict,
+    const objectRegistry& db
+)
 :
-    saturationModel(),
+    saturationModel(db),
     C_(dict.lookup("C<8>"))
 {}
 

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/polynomial/polynomial.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -77,7 +77,7 @@ public:
     // Constructors
 
         //- Construct from a dictionary
-        polynomial(const dictionary& dict);
+        polynomial(const dictionary& dict, const objectRegistry& db);
 
 
     //- Destructor

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/saturationModel/newSaturationModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,7 +29,8 @@ License
 
 Foam::autoPtr<Foam::saturationModel> Foam::saturationModel::New
 (
-    const dictionary& dict
+    const dictionary& dict,
+    const objectRegistry& db
 )
 {
     word saturationModelType(dict.lookup("type"));
@@ -50,7 +51,7 @@ Foam::autoPtr<Foam::saturationModel> Foam::saturationModel::New
             << exit(FatalError);
     }
 
-    return cstrIter()(dict);
+    return cstrIter()(dict, db);
 }
 
 

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/saturationModel/saturationModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -36,7 +36,17 @@ namespace Foam
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::saturationModel::saturationModel()
+Foam::saturationModel::saturationModel(const objectRegistry& db)
+:
+    IOdictionary
+    (
+        IOobject
+        (
+            "saturationModel",
+            db.time().constant(),
+            db
+        )
+    )
 {}
 
 

applications/solvers/multiphase/reactingEulerFoam/interfacialCompositionModels/saturationModels/saturationModel/saturationModel.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,6 +35,7 @@ SourceFiles
 #ifndef saturationModel_H
 #define saturationModel_H
 
+#include "IOdictionary.H"
 #include "volFields.H"
 #include "dictionary.H"
 #include "runTimeSelectionTables.H"
@@ -49,6 +50,8 @@ namespace Foam
 \*---------------------------------------------------------------------------*/
 
 class saturationModel
+:
+    public IOdictionary
 {
     // Private Member Functions
 
@@ -72,22 +75,26 @@ public:
         saturationModel,
         dictionary,
         (
-            const dictionary& dict
+            const dictionary& dict, const objectRegistry& db
         ),
-        (dict)
+        (dict, db)
     );
 
 
     // Constructors
 
         //- Construct null
-        saturationModel();
+        saturationModel(const objectRegistry& db);
 
 
     // Selectors
 
         //- Select null constructed
-        static autoPtr<saturationModel> New(const dictionary& dict);
+        static autoPtr<saturationModel> New
+        (
+            const dictionary& dict,
+            const objectRegistry& db
+        );
 
 
     //- Destructor

applications/solvers/multiphase/reactingEulerFoam/interfacialModels/heatTransferModels/heatTransferModel/heatTransferModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -53,7 +53,9 @@ Foam::heatTransferModel::heatTransferModel
         dict.lookupOrDefault<scalar>
         (
             "residualAlpha",
-            pair_.dispersed().residualAlpha().value()
+            pair_.ordered()
+          ? pair_.dispersed().residualAlpha().value()
+          : pair_.phase1().residualAlpha().value()
         )
     )
 {}

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -73,71 +73,50 @@ HeatAndMassTransferPhaseSystem
             continue;
         }
 
-        // Initialy assume no mass transfer
+        if
-
-        dmdt_.insert
         (
-            pair,
+            heatTransferModels_.found(pair)
-            new volScalarField
+         && heatTransferModels_[pair][pair.first()].valid()
-            (
+         && heatTransferModels_[pair][pair.second()].valid()
-                IOobject
+        )
-                (
+        {
-                    IOobject::groupName("dmdt", pair.name()),
-                    this->mesh().time().timeName(),
-                    this->mesh(),
-                    IOobject::NO_READ,
-                    IOobject::AUTO_WRITE
-                ),
-                this->mesh(),
-                dimensionedScalar("zero", dimDensity/dimTime, 0)
-            )
-        );
 
-        dmdtExplicit_.insert
+            volScalarField H1(heatTransferModels_[pair][pair.first()]->K());
-        (
+            volScalarField H2(heatTransferModels_[pair][pair.second()]->K());
-            pair,
-            new volScalarField
-            (
-                IOobject
-                (
-                    IOobject::groupName("dmdtExplicit", pair.name()),
-                    this->mesh().time().timeName(),
-                    this->mesh()
-                ),
-                this->mesh(),
-                dimensionedScalar("zero", dimDensity/dimTime, 0)
-            )
-        );
 
-        volScalarField H1(heatTransferModels_[pair][pair.first()]->K());
+            Tf_.insert
-        volScalarField H2(heatTransferModels_[pair][pair.second()]->K());
-
-        Tf_.insert
-        (
-            pair,
-            new volScalarField
             (
-                IOobject
+                pair,
+                new volScalarField
                 (
-                    IOobject::groupName("Tf", pair.name()),
+                    IOobject
-                    this->mesh().time().timeName(),
+                    (
-                    this->mesh(),
+                        IOobject::groupName("Tf", pair.name()),
-                    IOobject::NO_READ,
+                        this->mesh().time().timeName(),
-                    IOobject::AUTO_WRITE
+                        this->mesh(),
-                ),
+                        IOobject::NO_READ,
-                (
+                        IOobject::AUTO_WRITE
-                    H1*pair.phase1().thermo().T()
+                    ),
-                  + H2*pair.phase2().thermo().T()
+                    (
+                        H1*pair.phase1().thermo().T()
+                      + H2*pair.phase2().thermo().T()
+                    )
+                   /max
+                    (
+                        H1 + H2,
+                        dimensionedScalar
+                        (
+                            "small",
+                            heatTransferModel::dimK,
+                            SMALL
+                        )
+                    ),
+                    zeroGradientFvPatchScalarField::typeName
                 )
-               /max
+            );
-                (
+
-                    H1 + H2,
+            Tf_[pair]->correctBoundaryConditions();
-                    dimensionedScalar("small", heatTransferModel::dimK, SMALL)
+        }
-                ),
-                zeroGradientFvPatchScalarField::typeName
-            )
-        );
-        Tf_[pair]->correctBoundaryConditions();
     }
 }
 
@@ -169,9 +148,20 @@ Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::dmdt
     const phasePairKey& key
 ) const
 {
-    const scalar dmdtSign(Pair<word>::compare(dmdt_.find(key).key(), key));
+    return tmp<Foam::volScalarField>
-
+    (
-    return dmdtSign**dmdt_[key];
+        new volScalarField
+        (
+            IOobject
+            (
+                "dmdt",
+                this->mesh_.time().timeName(),
+                this->mesh_
+            ),
+            this->mesh_,
+            dimensionedScalar("zero", dimDensity/dimTime, 0)
+        )
+    );
 }
 
 
@@ -182,7 +172,7 @@ Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::dmdt
     const Foam::phaseModel& phase
 ) const
 {
-    tmp<volScalarField> tdmdt
+    tmp<volScalarField> tDmdt
     (
         new volScalarField
         (
@@ -197,35 +187,7 @@ Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::dmdt
         )
     );
 
-    forAllConstIter
+    return tDmdt;
-    (
-        phaseSystem::phasePairTable,
-        this->phasePairs_,
-        phasePairIter
-    )
-    {
-        const phasePair& pair(phasePairIter());
-
-        if (pair.ordered())
-        {
-            continue;
-        }
-
-        const phaseModel* phase1 = &pair.phase1();
-        const phaseModel* phase2 = &pair.phase2();
-
-        forAllConstIter(phasePair, pair, iter)
-        {
-            if (phase1 == &phase)
-            {
-                tdmdt.ref() += this->dmdt(pair);
-            }
-
-            Swap(phase1, phase2);
-        }
-    }
-
-    return tdmdt;
 }
 
 
@@ -236,34 +198,6 @@ Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::momentumTransfer() const
     autoPtr<phaseSystem::momentumTransferTable>
         eqnsPtr(BasePhaseSystem::momentumTransfer());
 
-    phaseSystem::momentumTransferTable& eqns = eqnsPtr();
-
-    // Source term due to mass trasfer
-    forAllConstIter
-    (
-        phaseSystem::phasePairTable,
-        this->phasePairs_,
-        phasePairIter
-    )
-    {
-        const phasePair& pair(phasePairIter());
-
-        if (pair.ordered())
-        {
-            continue;
-        }
-
-        const volVectorField& U1(pair.phase1().U());
-        const volVectorField& U2(pair.phase2().U());
-
-        const volScalarField dmdt(this->dmdt(pair));
-        const volScalarField dmdt21(posPart(dmdt));
-        const volScalarField dmdt12(negPart(dmdt));
-
-        *eqns[pair.phase1().name()] += dmdt21*U2 - fvm::Sp(dmdt21, U1);
-        *eqns[pair.phase2().name()] -= dmdt12*U1 - fvm::Sp(dmdt12, U2);
-    }
-
     return eqnsPtr;
 }
 
@@ -298,52 +232,98 @@ Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::heatTransfer() const
         heatTransferModelIter
     )
     {
-        const phasePair& pair
+        if
         (
-            this->phasePairs_[heatTransferModelIter.key()]
+            heatTransferModels_.found(heatTransferModelIter.key())
-        );
+        )
-
-        const phaseModel* phase = &pair.phase1();
-        const phaseModel* otherPhase = &pair.phase2();
-
-        const volScalarField& Tf(*Tf_[pair]);
-
-        const volScalarField K1
-        (
-            heatTransferModelIter()[pair.first()]->K()
-        );
-        const volScalarField K2
-        (
-            heatTransferModelIter()[pair.second()]->K()
-        );
-        const volScalarField KEff
-        (
-            K1*K2
-           /max
-            (
-                K1 + K2,
-                dimensionedScalar("small", heatTransferModel::dimK, SMALL)
-            )
-        );
-
-        const volScalarField* K = &K1;
-        const volScalarField* otherK = &K2;
-
-        forAllConstIter(phasePair, pair, iter)
         {
-            const volScalarField& he(phase->thermo().he());
+            const phasePair& pair
-            volScalarField Cpv(phase->thermo().Cpv());
+            (
+                this->phasePairs_[heatTransferModelIter.key()]
+            );
 
-            *eqns[phase->name()] +=
+            const phaseModel* phase = &pair.phase1();
-                (*K)*(Tf - phase->thermo().T())
+            const phaseModel* otherPhase = &pair.phase2();
-              + KEff/Cpv*he - fvm::Sp(KEff/Cpv, he);
 
-            Swap(phase, otherPhase);
+            const volScalarField& Tf(*Tf_[pair]);
-            Swap(K, otherK);
+
+            const volScalarField K1
+            (
+                heatTransferModelIter()[pair.first()]->K()
+            );
+            const volScalarField K2
+            (
+                heatTransferModelIter()[pair.second()]->K()
+            );
+            const volScalarField KEff
+            (
+                K1*K2
+               /max
+                (
+                    K1 + K2,
+                    dimensionedScalar("small", heatTransferModel::dimK, SMALL)
+                )
+            );
+
+            const volScalarField* K = &K1;
+            const volScalarField* otherK = &K2;
+
+            forAllConstIter(phasePair, pair, iter)
+            {
+                const volScalarField& he(phase->thermo().he());
+                volScalarField Cpv(phase->thermo().Cpv());
+
+                *eqns[phase->name()] +=
+                    (*K)*(Tf - phase->thermo().T())
+                  + KEff/Cpv*he - fvm::Sp(KEff/Cpv, he);
+
+                Swap(phase, otherPhase);
+                Swap(K, otherK);
+            }
         }
     }
 
-    // Source term due to mass transfer
+    return eqnsPtr;
+}
+
+template<class BasePhaseSystem>
+Foam::autoPtr<Foam::phaseSystem::massTransferTable>
+Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::massTransfer() const
+{
+    // Create a mass transfer matrix for each species of each phase
+    autoPtr<phaseSystem::massTransferTable> eqnsPtr
+    (
+        new phaseSystem::massTransferTable()
+    );
+
+    phaseSystem::massTransferTable& eqns = eqnsPtr();
+
+    forAll(this->phaseModels_, phasei)
+    {
+        const phaseModel& phase = this->phaseModels_[phasei];
+
+        const PtrList<volScalarField>& Yi = phase.Y();
+
+        forAll(Yi, i)
+        {
+            eqns.insert
+            (
+                Yi[i].name(),
+                new fvScalarMatrix(Yi[i], dimMass/dimTime)
+            );
+        }
+    }
+
+    return eqnsPtr;
+}
+
+
+template<class BasePhaseSystem>
+void Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::correctThermo()
+{
+
+    phaseSystem::correctThermo();
+
     forAllConstIter
     (
         phaseSystem::phasePairTable,
@@ -351,39 +331,49 @@ Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::heatTransfer() const
         phasePairIter
     )
     {
-        const phasePair& pair(phasePairIter());
+        if
-
+        (
-        if (pair.ordered())
+            this->heatTransferModels_.found(phasePairIter.key())
+        )
         {
-            continue;
+            const phasePair& pair(phasePairIter());
+
+            if (pair.ordered())
+            {
+                continue;
+            }
+
+            const phaseModel& phase1 = pair.phase1();
+            const phaseModel& phase2 = pair.phase2();
+
+            const volScalarField& T1(phase1.thermo().T());
+            const volScalarField& T2(phase2.thermo().T());
+
+            volScalarField& Tf(*this->Tf_[pair]);
+
+            volScalarField H1
+            (
+                this->heatTransferModels_[pair][pair.first()]->K()
+            );
+
+            volScalarField H2
+            (
+                this->heatTransferModels_[pair][pair.second()]->K()
+            );
+
+            // Limit the H[12] to avoid /0
+            H1.max(SMALL);
+            H2.max(SMALL);
+
+            Tf = (H1*T1 + H2*T2)/(H1 + H2);
+
+            Info<< "Tf." << pair.name()
+                << ": min = " << min(Tf.primitiveField())
+                << ", mean = " << average(Tf.primitiveField())
+                << ", max = " << max(Tf.primitiveField())
+                << endl;
         }
-
-        const phaseModel& phase1 = pair.phase1();
-        const phaseModel& phase2 = pair.phase2();
-
-        const volScalarField& he1(phase1.thermo().he());
-        const volScalarField& he2(phase2.thermo().he());
-
-        const volScalarField& K1(phase1.K());
-        const volScalarField& K2(phase2.K());
-
-        const volScalarField dmdt(this->dmdt(pair));
-        const volScalarField dmdt21(posPart(dmdt));
-        const volScalarField dmdt12(negPart(dmdt));
-        const volScalarField& Tf(*Tf_[pair]);
-
-        *eqns[phase1.name()] +=
-            dmdt21*(phase1.thermo().he(phase1.thermo().p(), Tf))
-          - fvm::Sp(dmdt21, he1)
-          + dmdt21*(K2 - K1);
-
-        *eqns[phase2.name()] -=
-            dmdt12*(phase2.thermo().he(phase2.thermo().p(), Tf))
-          - fvm::Sp(dmdt12, he2)
-          + dmdt12*(K1 - K2);
     }
-
-    return eqnsPtr;
 }
 
 

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -26,7 +26,8 @@ Class
 
 Description
     Base class to support interfacial heat and mass transfer between a number
-    of phases.
+    of phases. Can be used on its own to simulate flows without mass transfer
+    with two-resistance interfacial heat trasnfer.
 
 SourceFiles
     HeatAndMassTransferPhaseSystem.C
@@ -86,14 +87,6 @@ protected:
 
     // Protected data
 
-        //- Mass transfer rate
-        HashPtrTable<volScalarField, phasePairKey, phasePairKey::hash>
-            dmdt_;
-
-        //- Explicit part of the mass transfer rate
-        HashPtrTable<volScalarField, phasePairKey, phasePairKey::hash>
-            dmdtExplicit_;
-
         //- Interface temperatures
         HashPtrTable<volScalarField, phasePairKey, phasePairKey::hash> Tf_;
 
@@ -138,11 +131,10 @@ public:
             heatTransfer() const;
 
         //- Return the mass transfer matrices
-        virtual autoPtr<phaseSystem::massTransferTable>
+        virtual autoPtr<phaseSystem::massTransferTable> massTransfer() const;
-            massTransfer() const = 0;
 
         //- Correct the thermodynamics
-        virtual void correctThermo() = 0;
+        virtual void correctThermo();
 
         //- Read base phaseProperties dictionary
         virtual bool read();

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -106,8 +106,7 @@ public:
             heatTransfer() const;
 
         //- Return the mass transfer matrices
-        virtual autoPtr<phaseSystem::massTransferTable>
+        virtual autoPtr<phaseSystem::massTransferTable> massTransfer() const;
-            massTransfer() const;
 
         //- Read base phaseProperties dictionary
         virtual bool read();

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -36,13 +36,63 @@ InterfaceCompositionPhaseChangePhaseSystem
     const fvMesh& mesh
 )
 :
-    HeatAndMassTransferPhaseSystem<BasePhaseSystem>(mesh)
+    BasePhaseSystem(mesh)
 {
     this->generatePairsAndSubModels
     (
         "interfaceComposition",
         interfaceCompositionModels_
     );
+
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        const phasePair& pair(phasePairIter());
+
+        if (pair.ordered())
+        {
+            continue;
+        }
+
+        // Initially assume no mass transfer
+        iDmdt_.insert
+        (
+            pair,
+            new volScalarField
+            (
+                IOobject
+                (
+                    IOobject::groupName("iDmdt", pair.name()),
+                    this->mesh().time().timeName(),
+                    this->mesh(),
+                    IOobject::READ_IF_PRESENT,
+                    IOobject::AUTO_WRITE
+                ),
+                this->mesh(),
+                dimensionedScalar("zero", dimDensity/dimTime, 0)
+            )
+        );
+
+        iDmdtExplicit_.insert
+        (
+            pair,
+            new volScalarField
+            (
+                IOobject
+                (
+                    IOobject::groupName("iDmdtExplicit", pair.name()),
+                    this->mesh().time().timeName(),
+                    this->mesh()
+                ),
+                this->mesh(),
+                dimensionedScalar("zero", dimDensity/dimTime, 0)
+            )
+        );
+    }
 }
 
 
@@ -56,35 +106,242 @@ Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
 
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
+template<class BasePhaseSystem>
+Foam::tmp<Foam::volScalarField>
+Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::iDmdt
+(
+    const phasePairKey& key
+) const
+{
+    const scalar dmdtSign(Pair<word>::compare(iDmdt_.find(key).key(), key));
+
+    return dmdtSign**iDmdt_[key];
+}
+
+
+template<class BasePhaseSystem>
+Foam::tmp<Foam::volScalarField>
+Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::iDmdt
+(
+    const Foam::phaseModel& phase
+) const
+{
+    tmp<volScalarField> tiDmdt
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("iDmdt", phase.name()),
+                this->mesh_.time().timeName(),
+                this->mesh_
+            ),
+            this->mesh_,
+            dimensionedScalar("zero", dimDensity/dimTime, 0)
+        )
+    );
+
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        const phasePair& pair(phasePairIter());
+
+        if (pair.ordered())
+        {
+            continue;
+        }
+
+        const phaseModel* phase1 = &pair.phase1();
+        const phaseModel* phase2 = &pair.phase2();
+
+        forAllConstIter(phasePair, pair, iter)
+        {
+            if (phase1 == &phase)
+            {
+                tiDmdt.ref() += this->iDmdt
+                (
+                    phasePairKey(phase1->name(), phase2->name(),false)
+                );
+            }
+
+            Swap(phase1, phase2);
+        }
+    }
+
+    return tiDmdt;
+}
+
+
+template<class BasePhaseSystem>
+Foam::tmp<Foam::volScalarField>
+Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::dmdt
+(
+    const Foam::phaseModel& phase
+) const
+{
+    tmp<volScalarField> tDmdt
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("dmdt", phase.name()),
+                this->mesh_.time().timeName(),
+                this->mesh_
+            ),
+            this->mesh_,
+            dimensionedScalar("zero", dimDensity/dimTime, 0)
+        )
+    );
+
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        const phasePair& pair(phasePairIter());
+
+        if (pair.ordered())
+        {
+            continue;
+        }
+
+        const phaseModel* phase1 = &pair.phase1();
+        const phaseModel* phase2 = &pair.phase2();
+
+        forAllConstIter(phasePair, pair, iter)
+        {
+            if (phase1 == &phase)
+            {
+                tDmdt.ref() += this->dmdt
+                (
+                    phasePairKey(phase1->name(), phase2->name(),false)
+                );
+            }
+
+            Swap(phase1, phase2);
+        }
+    }
+    return tDmdt;
+}
+
+
+template<class BasePhaseSystem>
+Foam::autoPtr<Foam::phaseSystem::momentumTransferTable>
+Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
+momentumTransfer() const
+{
+    autoPtr<phaseSystem::momentumTransferTable>
+        eqnsPtr(BasePhaseSystem::momentumTransfer());
+
+    phaseSystem::momentumTransferTable& eqns = eqnsPtr();
+
+    // Source term due to mass transfer
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        const phasePair& pair(phasePairIter());
+
+        if (pair.ordered())
+        {
+            continue;
+        }
+
+        const volVectorField& U1(pair.phase1().U());
+        const volVectorField& U2(pair.phase2().U());
+
+        const volScalarField dmdt(this->iDmdt(pair));
+        const volScalarField dmdt21(posPart(dmdt));
+        const volScalarField dmdt12(negPart(dmdt));
+
+        *eqns[pair.phase1().name()] += dmdt21*U2 - fvm::Sp(dmdt21, U1);
+        *eqns[pair.phase2().name()] -= dmdt12*U1 - fvm::Sp(dmdt12, U2);
+    }
+
+    return eqnsPtr;
+}
+
+
+template<class BasePhaseSystem>
+Foam::autoPtr<Foam::phaseSystem::heatTransferTable>
+Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
+heatTransfer() const
+{
+    autoPtr<phaseSystem::heatTransferTable> eqnsPtr =
+        BasePhaseSystem::heatTransfer();
+
+    phaseSystem::heatTransferTable& eqns = eqnsPtr();
+
+    // Source term due to mass transfer
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        if
+        (
+            this->heatTransferModels_.found(phasePairIter.key())
+        )
+        {
+            const phasePair& pair(phasePairIter());
+
+            if (pair.ordered())
+            {
+                continue;
+            }
+
+            const phaseModel& phase1 = pair.phase1();
+            const phaseModel& phase2 = pair.phase2();
+
+            const volScalarField& he1(phase1.thermo().he());
+            const volScalarField& he2(phase2.thermo().he());
+
+            const volScalarField& K1(phase1.K());
+            const volScalarField& K2(phase2.K());
+
+            const volScalarField dmdt(this->dmdt(pair));
+            const volScalarField dmdt21(posPart(dmdt));
+            const volScalarField dmdt12(negPart(dmdt));
+            const volScalarField& Tf(*this->Tf_[pair]);
+
+            *eqns[phase1.name()] +=
+                dmdt21*(phase1.thermo().he(phase1.thermo().p(), Tf))
+              - fvm::Sp(dmdt21, he1)
+              + dmdt21*(K2 - K1);
+
+            *eqns[phase2.name()] -=
+                dmdt12*(phase2.thermo().he(phase2.thermo().p(), Tf))
+              - fvm::Sp(dmdt12, he2)
+              + dmdt12*(K1 - K2);
+        }
+    }
+
+    return eqnsPtr;
+}
+
+
 template<class BasePhaseSystem>
 Foam::autoPtr<Foam::phaseSystem::massTransferTable>
 Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
 massTransfer() const
 {
-    // Create a mass transfer matrix for each species of each phase
+    autoPtr<phaseSystem::massTransferTable> eqnsPtr =
-    autoPtr<phaseSystem::massTransferTable> eqnsPtr
+        BasePhaseSystem::massTransfer();
-    (
-        new phaseSystem::massTransferTable()
-    );
 
     phaseSystem::massTransferTable& eqns = eqnsPtr();
 
-    forAll(this->phaseModels_, phasei)
-    {
-        const phaseModel& phase = this->phaseModels_[phasei];
-
-        const PtrList<volScalarField>& Yi = phase.Y();
-
-        forAll(Yi, i)
-        {
-            eqns.insert
-            (
-                Yi[i].name(),
-                new fvScalarMatrix(Yi[i], dimMass/dimTime)
-            );
-        }
-    }
-
     // Reset the interfacial mass flow rates
     forAllConstIter
     (
@@ -100,10 +357,10 @@ massTransfer() const
             continue;
         }
 
-        *this->dmdt_[pair] =
+        *this->iDmdt_[pair] =
-            *this->dmdtExplicit_[pair];
+            *this->iDmdtExplicit_[pair];
 
-        *this->dmdtExplicit_[pair] =
+        *this->iDmdtExplicit_[pair] =
             dimensionedScalar("zero", dimDensity/dimTime, 0);
     }
 
@@ -130,10 +387,10 @@ massTransfer() const
 
         const volScalarField& Tf(*this->Tf_[key]);
 
-        volScalarField& dmdtExplicit(*this->dmdtExplicit_[key]);
+        volScalarField& iDmdtExplicit_(*this->iDmdtExplicit_[key]);
-        volScalarField& dmdt(*this->dmdt_[key]);
+        volScalarField& iDmdt_(*this->iDmdt_[key]);
 
-        scalar dmdtSign(Pair<word>::compare(this->dmdt_.find(key).key(), key));
+        scalar dmdtSign(Pair<word>::compare(this->iDmdt_.find(key).key(), key));
 
         const volScalarField K
         (
@@ -175,8 +432,8 @@ massTransfer() const
               - fvm::Sp(phase.rho()*KD, eqns[name]->psi());
 
             // Sum the mass transfer rate
-            dmdtExplicit += dmdtSign*phase.rho()*KD*Yf;
+            iDmdtExplicit_ += dmdtSign*phase.rho()*KD*Yf;
-            dmdt -= dmdtSign*phase.rho()*KD*eqns[name]->psi();
+            iDmdt_ -= dmdtSign*phase.rho()*KD*eqns[name]->psi();
 
             // Explicit transport out of the other phase
             if (eqns.found(otherName))
@@ -195,7 +452,7 @@ template<class BasePhaseSystem>
 void Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
 correctThermo()
 {
-    BasePhaseSystem::correctThermo();
+    phaseSystem::correctThermo();
 
     // This loop solves for the interface temperatures, Tf, and updates the
     // interface composition models.

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -40,7 +40,7 @@ SourceFiles
 #ifndef InterfaceCompositionPhaseChangePhaseSystem_H
 #define InterfaceCompositionPhaseChangePhaseSystem_H
 
-#include "HeatAndMassTransferPhaseSystem.H"
+#include "phaseSystem.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -56,7 +56,7 @@ class interfaceCompositionModel;
 template<class BasePhaseSystem>
 class InterfaceCompositionPhaseChangePhaseSystem
 :
-    public HeatAndMassTransferPhaseSystem<BasePhaseSystem>
+    public BasePhaseSystem
 {
 protected:
 
@@ -77,6 +77,14 @@ protected:
             //- Interface composition models
             interfaceCompositionModelTable interfaceCompositionModels_;
 
+            //- Interfacial Mass transfer rate
+            HashPtrTable<volScalarField, phasePairKey, phasePairKey::hash>
+                iDmdt_;
+
+            //- Explicit part of the mass transfer rate
+            HashPtrTable<volScalarField, phasePairKey, phasePairKey::hash>
+                iDmdtExplicit_;
+
 
 public:
 
@@ -92,6 +100,28 @@ public:
 
     // Member Functions
 
+        //- Return the interfacial mass flow rate
+        virtual tmp<volScalarField> iDmdt(const phasePairKey& key) const;
+
+        //- Return the total interfacial mass transfer rate for phase
+        virtual tmp<volScalarField> iDmdt(const phaseModel& phase) const;
+
+        //- Return the interfacial mass flow rate
+        virtual tmp<volScalarField> dmdt(const phasePairKey& key) const
+        {
+            return BasePhaseSystem::dmdt(key) + this->iDmdt(key);
+        };
+
+        //- Return the total interfacial mass transfer rate for phase
+        virtual tmp<volScalarField> dmdt(const phaseModel& phase) const;
+
+        //- Return the momentum transfer matrices
+        virtual autoPtr<phaseSystem::momentumTransferTable>
+            momentumTransfer() const;
+
+        //- Return the heat transfer matrices
+        virtual autoPtr<phaseSystem::heatTransferTable> heatTransfer() const;
+
         //- Return the mass transfer matrices
         virtual autoPtr<phaseSystem::massTransferTable> massTransfer() const;
 

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C

@@ -26,6 +26,7 @@ License
 #include "ThermalPhaseChangePhaseSystem.H"
 #include "alphatPhaseChangeWallFunctionFvPatchScalarField.H"
 #include "fvcVolumeIntegrate.H"
+#include "fvmSup.H"
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
@@ -36,9 +37,12 @@ ThermalPhaseChangePhaseSystem
     const fvMesh& mesh
 )
 :
-    HeatAndMassTransferPhaseSystem<BasePhaseSystem>(mesh),
+    BasePhaseSystem(mesh),
-    volatile_(this->lookup("volatile")),
+    volatile_(this->template lookupOrDefault<word>("volatile","none")),
-    saturationModel_(saturationModel::New(this->subDict("saturationModel"))),
+    saturationModel_
+    (
+        saturationModel::New(this->subDict("saturationModel"), mesh)
+    ),
     massTransfer_(this->lookup("massTransfer"))
 {
 
@@ -74,6 +78,44 @@ ThermalPhaseChangePhaseSystem
                 dimensionedScalar("zero", dimDensity/dimTime, 0)
             )
         );
+
+        // Initially assume no mass transfer
+        wDmdt_.insert
+        (
+            pair,
+            new volScalarField
+            (
+                IOobject
+                (
+                    IOobject::groupName("wDmdt", pair.name()),
+                    this->mesh().time().timeName(),
+                    this->mesh(),
+                    IOobject::READ_IF_PRESENT,
+                    IOobject::AUTO_WRITE
+                ),
+                this->mesh(),
+                dimensionedScalar("zero", dimDensity/dimTime, 0)
+            )
+        );
+
+        // Initially assume no mass transfer
+        wMDotL_.insert
+        (
+            pair,
+            new volScalarField
+            (
+                IOobject
+                (
+                    IOobject::groupName("wMDotL", pair.name()),
+                    this->mesh().time().timeName(),
+                    this->mesh(),
+                    IOobject::READ_IF_PRESENT,
+                    IOobject::AUTO_WRITE
+                ),
+                this->mesh(),
+                dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0)
+            )
+        );
     }
 }
 
@@ -96,167 +138,6 @@ Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::saturation() const
 }
 
 
-template<class BasePhaseSystem>
-Foam::autoPtr<Foam::phaseSystem::heatTransferTable>
-Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::heatTransfer() const
-{
-    typedef compressible::alphatPhaseChangeWallFunctionFvPatchScalarField
-        alphatPhaseChangeWallFunction;
-
-    autoPtr<phaseSystem::heatTransferTable> eqnsPtr =
-        Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::heatTransfer();
-
-    phaseSystem::heatTransferTable& eqns = eqnsPtr();
-
-    // Accumulate mDotL contributions from boundaries
-    forAllConstIter
-    (
-        phaseSystem::phasePairTable,
-        this->phasePairs_,
-        phasePairIter
-    )
-    {
-        const phasePair& pair(phasePairIter());
-
-        if (pair.ordered())
-        {
-            continue;
-        }
-
-        const phaseModel& phase = pair.phase1();
-        const phaseModel& otherPhase = pair.phase2();
-
-        volScalarField mDotL
-        (
-            IOobject
-            (
-                "mDotL",
-                phase.mesh().time().timeName(),
-                phase.mesh(),
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            phase.mesh(),
-            dimensionedScalar("",dimensionSet(1,-1,-3,0,0),0.0)
-        );
-
-        if
-        (
-            otherPhase.mesh().foundObject<volScalarField>
-            (
-                "alphat." +  otherPhase.name()
-            )
-        )
-        {
-            const volScalarField& alphat =
-                otherPhase.mesh().lookupObject<volScalarField>
-                (
-                    "alphat." +  otherPhase.name()
-                );
-
-            const fvPatchList& patches = this->mesh().boundary();
-            forAll(patches, patchi)
-            {
-                const fvPatch& currPatch = patches[patchi];
-
-                if
-                (
-                    isA<alphatPhaseChangeWallFunction>
-                    (
-                        alphat.boundaryField()[patchi]
-                    )
-                )
-                {
-                    const scalarField& patchMDotL =
-                        refCast<const alphatPhaseChangeWallFunction>
-                        (
-                            alphat.boundaryField()[patchi]
-                        ).mDotL();
-
-                    forAll(patchMDotL,facei)
-                    {
-                        label faceCelli = currPatch.faceCells()[facei];
-                        mDotL[faceCelli] = patchMDotL[facei];
-                    }
-                }
-            }
-        }
-
-        *eqns[otherPhase.name()] -= mDotL;
-
-    }
-
-    return eqnsPtr;
-}
-
-
-template<class BasePhaseSystem>
-Foam::autoPtr<Foam::phaseSystem::massTransferTable>
-Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::massTransfer() const
-{
-    // Create a mass transfer matrix for each species of each phase
-    autoPtr<phaseSystem::massTransferTable> eqnsPtr
-    (
-        new phaseSystem::massTransferTable()
-    );
-
-    phaseSystem::massTransferTable& eqns = eqnsPtr();
-
-    forAll(this->phaseModels_, phasei)
-    {
-        const phaseModel& phase = this->phaseModels_[phasei];
-
-        const PtrList<volScalarField>& Yi = phase.Y();
-
-        forAll(Yi, i)
-        {
-            eqns.insert
-            (
-                Yi[i].name(),
-                new fvScalarMatrix(Yi[i], dimMass/dimTime)
-            );
-        }
-    }
-
-    forAllConstIter
-    (
-        phaseSystem::phasePairTable,
-        this->phasePairs_,
-        phasePairIter
-    )
-    {
-        const phasePair& pair(phasePairIter());
-
-        if (pair.ordered())
-        {
-            continue;
-        }
-        const phaseModel& phase = pair.phase1();
-        const phaseModel& otherPhase = pair.phase2();
-
-        const word name
-        (
-            IOobject::groupName(volatile_, phase.name())
-        );
-
-        const word otherName
-        (
-            IOobject::groupName(volatile_, otherPhase.name())
-        );
-
-        const volScalarField dmdt(this->dmdt(pair));
-        const volScalarField dmdt12(posPart(dmdt));
-        const volScalarField dmdt21(negPart(dmdt));
-
-        *eqns[name] += fvm::Sp(dmdt21, eqns[name]->psi()) - dmdt21;
-        *eqns[otherName] += dmdt12 - fvm::Sp(dmdt12, eqns[otherName]->psi());
-    }
-
-    return eqnsPtr;
-}
-
-
 template<class BasePhaseSystem>
 Foam::tmp<Foam::volScalarField>
 Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::iDmdt
@@ -265,7 +146,6 @@ Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::iDmdt
 ) const
 {
     const scalar dmdtSign(Pair<word>::compare(iDmdt_.find(key).key(), key));
-
     return dmdtSign**iDmdt_[key];
 }
 
@@ -313,7 +193,10 @@ Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::iDmdt
         {
             if (phase1 == &phase)
             {
-                tiDmdt.ref() += this->iDmdt(pair);
+                tiDmdt.ref() += this->iDmdt
+                (
+                    phasePairKey(phase1->name(), phase2->name(),false)
+                );
             }
 
             Swap(phase1, phase2);
@@ -323,14 +206,40 @@ Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::iDmdt
     return tiDmdt;
 }
 
+template<class BasePhaseSystem>
+Foam::tmp<Foam::volScalarField>
+Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::wDmdt
+(
+    const phasePairKey& key
+) const
+{
+    const scalar dmdtSign(Pair<word>::compare(wDmdt_.find(key).key(), key));
+
+    return dmdtSign**wDmdt_[key];
+}
+
 
 template<class BasePhaseSystem>
-void Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::correctThermo()
+Foam::tmp<Foam::volScalarField>
+Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::wDmdt
+(
+    const Foam::phaseModel& phase
+) const
 {
-    typedef compressible::alphatPhaseChangeWallFunctionFvPatchScalarField
+    tmp<volScalarField> twDmdt
-        alphatPhaseChangeWallFunction;
+    (
-
+        new volScalarField
-    BasePhaseSystem::correctThermo();
+        (
+            IOobject
+            (
+                IOobject::groupName("wDmdt", phase.name()),
+                this->mesh_.time().timeName(),
+                this->mesh_
+            ),
+            this->mesh_,
+            dimensionedScalar("zero", dimDensity/dimTime, 0)
+        )
+    );
 
     forAllConstIter
     (
@@ -346,174 +255,437 @@ void Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::correctThermo()
             continue;
         }
 
-        const phaseModel& phase1 = pair.phase1();
+        const phaseModel* phase1 = &pair.phase1();
-        const phaseModel& phase2 = pair.phase2();
+        const phaseModel* phase2 = &pair.phase2();
 
-        Info<< phase1.name() << " min/max T "
+        forAllConstIter(phasePair, pair, iter)
-            << min(phase1.thermo().T()).value()
-            << " - "
-            << max(phase1.thermo().T()).value()
-            << endl;
-
-        Info<< phase2.name() << " min/max T "
-            << min(phase2.thermo().T()).value()
-            << " - "
-            << max(phase2.thermo().T()).value()
-            << endl;
-
-        const volScalarField& T1(phase1.thermo().T());
-        const volScalarField& T2(phase2.thermo().T());
-
-        const volScalarField& he1(phase1.thermo().he());
-        const volScalarField& he2(phase2.thermo().he());
-
-        volScalarField& dmdt(*this->dmdt_[pair]);
-        volScalarField& iDmdt(*this->iDmdt_[pair]);
-
-        volScalarField& Tf = *this->Tf_[pair];
-
-        volScalarField hef1(phase1.thermo().he(phase1.thermo().p(), Tf));
-        volScalarField hef2(phase2.thermo().he(phase2.thermo().p(), Tf));
-
-        volScalarField L
-        (
-            min
-            (
-                (pos0(iDmdt)*he2 + neg(iDmdt)*hef2)
-              - (neg(iDmdt)*he1 + pos0(iDmdt)*hef1),
-                0.3*mag(hef2 - hef1)
-            )
-        );
-
-        volScalarField iDmdtNew(iDmdt);
-
-        if (massTransfer_ )
         {
+            if (phase1 == &phase)
+            {
+                twDmdt.ref() += this->wDmdt
+                (
+                    phasePairKey(phase1->name(), phase2->name(),false)
+                );
+            }
+
+            Swap(phase1, phase2);
+        }
+    }
+
+    return twDmdt;
+}
+
+
+template<class BasePhaseSystem>
+Foam::tmp<Foam::volScalarField>
+Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::dmdt
+(
+    const Foam::phaseModel& phase
+) const
+{
+    tmp<volScalarField> tDmdt
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("Dmdt", phase.name()),
+                this->mesh_.time().timeName(),
+                this->mesh_
+            ),
+            this->mesh_,
+            dimensionedScalar("zero", dimDensity/dimTime, 0)
+        )
+    );
+
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        const phasePair& pair(phasePairIter());
+
+        if (pair.ordered())
+        {
+            continue;
+        }
+
+        const phaseModel* phase1 = &pair.phase1();
+        const phaseModel* phase2 = &pair.phase2();
+
+        forAllConstIter(phasePair, pair, iter)
+        {
+            if (phase1 == &phase)
+            {
+                tDmdt.ref() += this->dmdt
+                (
+                    phasePairKey(phase1->name(), phase2->name(),false)
+                );
+            }
+
+            Swap(phase1, phase2);
+        }
+    }
+
+    return tDmdt;
+}
+
+
+template<class BasePhaseSystem>
+Foam::autoPtr<Foam::phaseSystem::momentumTransferTable>
+Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::momentumTransfer() const
+{
+    autoPtr<phaseSystem::momentumTransferTable>
+        eqnsPtr(BasePhaseSystem::momentumTransfer());
+
+    phaseSystem::momentumTransferTable& eqns = eqnsPtr();
+
+    // Source term due to mass trasfer
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        const phasePair& pair(phasePairIter());
+
+        if (pair.ordered())
+        {
+            continue;
+        }
+
+        const volVectorField& U1(pair.phase1().U());
+        const volVectorField& U2(pair.phase2().U());
+
+        const volScalarField dmdt(this->iDmdt(pair) + this->wDmdt(pair));
+        const volScalarField dmdt21(posPart(dmdt));
+        const volScalarField dmdt12(negPart(dmdt));
+
+        *eqns[pair.phase1().name()] += dmdt21*U2 - fvm::Sp(dmdt21, U1);
+        *eqns[pair.phase2().name()] -= dmdt12*U1 - fvm::Sp(dmdt12, U2);
+    }
+
+    return eqnsPtr;
+}
+
+
+template<class BasePhaseSystem>
+Foam::autoPtr<Foam::phaseSystem::heatTransferTable>
+Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::heatTransfer() const
+{
+    autoPtr<phaseSystem::heatTransferTable> eqnsPtr =
+        BasePhaseSystem::heatTransfer();
+
+    phaseSystem::heatTransferTable& eqns = eqnsPtr();
+
+    // Source term due to mass transfer
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        if
+        (
+            this->heatTransferModels_.found(phasePairIter.key())
+        )
+        {
+            const phasePair& pair(phasePairIter());
+
+            if (pair.ordered())
+            {
+                continue;
+            }
+
+            const phaseModel& phase1 = pair.phase1();
+            const phaseModel& phase2 = pair.phase2();
+
+            const volScalarField& he1(phase1.thermo().he());
+            const volScalarField& he2(phase2.thermo().he());
+
+            const volScalarField& K1(phase1.K());
+            const volScalarField& K2(phase2.K());
+
+            const volScalarField dmdt(this->dmdt(pair));
+            const volScalarField dmdt21(posPart(dmdt));
+            const volScalarField dmdt12(negPart(dmdt));
+            const volScalarField& Tf(*this->Tf_[pair]);
+
+            *eqns[phase1.name()] +=
+                dmdt21*(phase1.thermo().he(phase1.thermo().p(), Tf))
+              - fvm::Sp(dmdt21, he1)
+              + dmdt21*(K2 - K1);
+
+            *eqns[phase2.name()] -=
+                dmdt12*(phase2.thermo().he(phase2.thermo().p(), Tf))
+              - fvm::Sp(dmdt12, he2)
+              + dmdt12*(K1 - K2);
+
+            *eqns[phase1.name()] += negPart(*this->wMDotL_[pair]);
+            *eqns[phase2.name()] -= posPart(*this->wMDotL_[pair]);
+        }
+    }
+
+    return eqnsPtr;
+}
+
+
+template<class BasePhaseSystem>
+Foam::autoPtr<Foam::phaseSystem::massTransferTable>
+Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::massTransfer() const
+{
+    autoPtr<phaseSystem::massTransferTable> eqnsPtr =
+        BasePhaseSystem::massTransfer();
+
+    phaseSystem::massTransferTable& eqns = eqnsPtr();
+
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        const phasePair& pair(phasePairIter());
+
+        if (pair.ordered())
+        {
+            continue;
+        }
+        const phaseModel& phase = pair.phase1();
+        const phaseModel& otherPhase = pair.phase2();
+
+        const PtrList<volScalarField>& Yi = phase.Y();
+
+        forAll(Yi, i)
+        {
+
+            if (Yi[i].member() != volatile_)
+            {
+                continue;
+            }
+
+            const word name
+            (
+                IOobject::groupName(volatile_, phase.name())
+            );
+
+            const word otherName
+            (
+                IOobject::groupName(volatile_, otherPhase.name())
+            );
+
+            const volScalarField dmdt(this->iDmdt(pair) + this->wDmdt(pair));
+            const volScalarField dmdt12(posPart(dmdt));
+            const volScalarField dmdt21(negPart(dmdt));
+
+            *eqns[name] += fvm::Sp(dmdt21, eqns[name]->psi()) - dmdt21;
+            *eqns[otherName] +=
+                dmdt12 - fvm::Sp(dmdt12, eqns[otherName]->psi());
+        }
+    }
+
+    return eqnsPtr;
+}
+
+
+template<class BasePhaseSystem>
+void Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::correctThermo()
+{
+    typedef compressible::alphatPhaseChangeWallFunctionFvPatchScalarField
+        alphatPhaseChangeWallFunction;
+
+    phaseSystem::correctThermo();
+
+    forAllConstIter
+    (
+        phaseSystem::phasePairTable,
+        this->phasePairs_,
+        phasePairIter
+    )
+    {
+        if
+        (
+            this->heatTransferModels_.found(phasePairIter.key())
+        )
+        {
+            const phasePair& pair(phasePairIter());
+
+            if (pair.ordered())
+            {
+                continue;
+            }
+
+            const phaseModel& phase1 = pair.phase1();
+            const phaseModel& phase2 = pair.phase2();
+
+            const volScalarField& T1(phase1.thermo().T());
+            const volScalarField& T2(phase2.thermo().T());
+
+            const volScalarField& he1(phase1.thermo().he());
+            const volScalarField& he2(phase2.thermo().he());
+
+            volScalarField& iDmdt(*this->iDmdt_[pair]);
+            volScalarField& Tf(*this->Tf_[pair]);
+
+            volScalarField hef1(phase1.thermo().he(phase1.thermo().p(), Tf));
+            volScalarField hef2(phase2.thermo().he(phase2.thermo().p(), Tf));
+
+            volScalarField L
+            (
+                (neg0(iDmdt)*hef2 + pos(iDmdt)*he2)
+              - (pos0(iDmdt)*hef1 + neg(iDmdt)*he1)
+            );
+
+            volScalarField iDmdtNew(iDmdt);
+
+            if (massTransfer_ )
+            {
+                volScalarField H1
+                (
+                    this->heatTransferModels_[pair][pair.first()]->K(0)
+                );
+
+                volScalarField H2
+                (
+                    this->heatTransferModels_[pair][pair.second()]->K(0)
+                );
+
+                Tf = saturationModel_->Tsat(phase1.thermo().p());
+
+                iDmdtNew =
+                    (H1*(Tf - T1) + H2*(Tf - T2))/L;
+            }
+            else
+            {
+                iDmdtNew == dimensionedScalar("0",iDmdt.dimensions(), 0);
+            }
+
             volScalarField H1
             (
-                this->heatTransferModels_[pair][pair.first()]->K(0)
+                this->heatTransferModels_[pair][pair.first()]->K()
             );
 
             volScalarField H2
             (
-                this->heatTransferModels_[pair][pair.second()]->K(0)
+                this->heatTransferModels_[pair][pair.second()]->K()
             );
 
-            Tf = saturationModel_->Tsat(phase1.thermo().p());
+            // Limit the H[12] to avoid /0
+            H1.max(SMALL);
+            H2.max(SMALL);
 
-            iDmdtNew =
+            Tf = (H1*T1 + H2*T2 + iDmdtNew*L)/(H1 + H2);
-                (H1*(Tf - T1) + H2*(Tf - T2))/L;
 
-        }
+            Info<< "Tf." << pair.name()
-        else
+                << ": min = " << min(Tf.primitiveField())
-        {
+                << ", mean = " << average(Tf.primitiveField())
-            iDmdtNew == dimensionedScalar("0",dmdt.dimensions(), 0);
+                << ", max = " << max(Tf.primitiveField())
-        }
-
-        volScalarField H1(this->heatTransferModels_[pair][pair.first()]->K());
-        volScalarField H2(this->heatTransferModels_[pair][pair.second()]->K());
-
-        // Limit the H[12] boundary field to avoid /0
-        const scalar HLimit = 1e-4;
-        H1.boundaryFieldRef() =
-            max(H1.boundaryField(), phase1.boundaryField()*HLimit);
-        H2.boundaryFieldRef() =
-            max(H2.boundaryField(), phase2.boundaryField()*HLimit);
-
-        Tf = (H1*T1 + H2*T2 + iDmdtNew*L)/(H1 + H2);
-
-        Info<< "Tf." << pair.name()
-            << ": min = " << min(Tf.primitiveField())
-            << ", mean = " << average(Tf.primitiveField())
-            << ", max = " << max(Tf.primitiveField())
-            << endl;
-
-        scalar iDmdtRelax(this->mesh().fieldRelaxationFactor("iDmdt"));
-        iDmdt = (1 - iDmdtRelax)*iDmdt + iDmdtRelax*iDmdtNew;
-
-        if (massTransfer_ )
-        {
-            Info<< "iDmdt." << pair.name()
-                << ": min = " << min(iDmdt.primitiveField())
-                << ", mean = " << average(iDmdt.primitiveField())
-                << ", max = " << max(iDmdt.primitiveField())
-                << ", integral = " << fvc::domainIntegrate(iDmdt).value()
                 << endl;
-        }
 
-        // Accumulate dmdt contributions from boundaries
+            scalar iDmdtRelax(this->mesh().fieldRelaxationFactor("iDmdt"));
-        volScalarField wDmdt
+            iDmdt = (1 - iDmdtRelax)*iDmdt + iDmdtRelax*iDmdtNew;
-        (
-            IOobject
-            (
-                IOobject::groupName("wDmdt", pair.name()),
-                this->mesh().time().timeName(),
-                this->mesh(),
-                IOobject::NO_READ,
-                IOobject::AUTO_WRITE,
-                false
-            ),
-            this->mesh(),
-            dimensionedScalar("zero", dimDensity/dimTime, 0)
-        );
 
-        if
+            if (massTransfer_ )
-        (
-            phase2.mesh().foundObject<volScalarField>
-            (
-                "alphat." +  phase2.name()
-            )
-        )
-        {
-            const volScalarField& alphat =
-                phase2.mesh().lookupObject<volScalarField>
-                (
-                    "alphat." +  phase2.name()
-                );
-
-            const fvPatchList& patches = this->mesh().boundary();
-            forAll(patches, patchi)
             {
-                const fvPatch& currPatch = patches[patchi];
+                Info<< "iDmdt." << pair.name()
+                    << ": min = " << min(iDmdt.primitiveField())
+                    << ", mean = " << average(iDmdt.primitiveField())
+                    << ", max = " << max(iDmdt.primitiveField())
+                    << ", integral = " << fvc::domainIntegrate(iDmdt).value()
+                    << endl;
+            }
 
+            volScalarField& wDmdt(*this->wDmdt_[pair]);
+            volScalarField& wMDotL(*this->wMDotL_[pair]);
+            wDmdt *= 0;
+            wMDotL *= 0;
+
+            bool wallBoilingActive = false;
+            const phaseModel* phase1Ptr = &phase1;
+            const phaseModel* phase2Ptr = &phase2;
+
+            forAllConstIter(phasePair, pair, iter)
+            {
                 if
                 (
-                    isA<alphatPhaseChangeWallFunction>
+                    phase1Ptr->mesh().foundObject<volScalarField>
                     (
-                        alphat.boundaryField()[patchi]
+                        "alphat." +  phase1Ptr->name()
                     )
                 )
                 {
-                    const scalarField& patchDmdt =
+                    const volScalarField& alphat =
-                        refCast<const alphatPhaseChangeWallFunction>
+                        phase1Ptr->mesh().lookupObject<volScalarField>
                         (
-                            alphat.boundaryField()[patchi]
+                            "alphat." +  phase1Ptr->name()
-                        ).dmdt();
+                        );
 
-                    forAll(patchDmdt,facei)
+                    const fvPatchList& patches = this->mesh().boundary();
+                    forAll(patches, patchi)
                     {
-                        label faceCelli = currPatch.faceCells()[facei];
+                        const fvPatch& currPatch = patches[patchi];
-                        wDmdt[faceCelli] += patchDmdt[facei];
+
+                        if
+                        (
+                            isA<alphatPhaseChangeWallFunction>
+                            (
+                                alphat.boundaryField()[patchi]
+                            )
+                        )
+                        {
+                            const alphatPhaseChangeWallFunction& PCpatch =
+                                refCast<const alphatPhaseChangeWallFunction>
+                                (
+                                    alphat.boundaryField()[patchi]
+                                );
+                            phasePairKey key
+                            (
+                                phase1Ptr->name(), phase2Ptr->name(), false
+                            );
+                            if(PCpatch.activePhasePair(key))
+                            {
+                                wallBoilingActive = true;
+
+                                const scalarField& patchDmdt =
+                                    PCpatch.dmdt(key);
+                                const scalarField& patchMDotL =
+                                    PCpatch.mDotL(key);
+
+                                const scalar sign
+                                (
+                                    Pair<word>::compare(pair, key)
+                                );
+
+                                forAll(patchDmdt,facei)
+                                {
+                                    const label faceCelli =
+                                        currPatch.faceCells()[facei];
+                                    wDmdt[faceCelli] -= sign*patchDmdt[facei];
+                                    wMDotL[faceCelli] -= sign*patchMDotL[facei];
+                                }
+                            }
+                        }
                     }
                 }
+
+                Swap(phase1Ptr, phase2Ptr);
             }
 
-            Info<< "wDmdt." << pair.name()
+            if (wallBoilingActive)
-                << ": min = " << min(wDmdt.primitiveField())
+            {
-                << ", mean = " << average(wDmdt.primitiveField())
+                Info<< "wDmdt." << pair.name()
-                << ", max = " << max(wDmdt.primitiveField())
+                    << ": min = " << min(wDmdt.primitiveField())
-                << ", integral = " << fvc::domainIntegrate(wDmdt).value()
+                    << ", mean = " << average(wDmdt.primitiveField())
-                << endl;
+                    << ", max = " << max(wDmdt.primitiveField())
+                    << ", integral = " << fvc::domainIntegrate(wDmdt).value()
+                    << endl;
+            }
         }
-
-        dmdt = wDmdt + iDmdt;
-
-        Info<< "dmdt." << pair.name()
-            << ": min = " << min(dmdt.primitiveField())
-            << ", mean = " << average(dmdt.primitiveField())
-            << ", max = " << max(dmdt.primitiveField())
-            << ", integral = " << fvc::domainIntegrate(dmdt).value()
-            << endl;
     }
 }
 

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,10 +29,6 @@ Description
     phases according the interfacial temperature approximated by the saturation
     temperature.
 
-    Currently only a single specified specie is considered volatile and changes
-    phase, all other species are considered nonvolatile and do not
-    affect the mass-transfer.
-
 SourceFiles
     ThermalPhaseChangePhaseSystem.C
 
@@ -41,7 +37,7 @@ SourceFiles
 #ifndef ThermalPhaseChangePhaseSystem_H
 #define ThermalPhaseChangePhaseSystem_H
 
-#include "HeatAndMassTransferPhaseSystem.H"
+#include "phaseSystem.H"
 #include "saturationModel.H"
 #include "Switch.H"
 
@@ -57,7 +53,7 @@ namespace Foam
 template<class BasePhaseSystem>
 class ThermalPhaseChangePhaseSystem
 :
-    public HeatAndMassTransferPhaseSystem<BasePhaseSystem>
+    public BasePhaseSystem
 {
 
 protected:
@@ -77,6 +73,15 @@ protected:
         HashPtrTable<volScalarField, phasePairKey, phasePairKey::hash>
             iDmdt_;
 
+        //- Boundary Mass transfer rate
+        HashPtrTable<volScalarField, phasePairKey, phasePairKey::hash>
+            wDmdt_;
+
+        //- Boundary thermal energy transfer rate
+        HashPtrTable<volScalarField, phasePairKey, phasePairKey::hash>
+            wMDotL_;
+
+
 
 public:
 
@@ -95,18 +100,38 @@ public:
         //- Return the saturationModel
         const saturationModel& saturation() const;
 
-        //- Return the heat transfer matrices
-        virtual autoPtr<phaseSystem::heatTransferTable> heatTransfer() const;
-
-        //- Return the mass transfer matrices
-        virtual autoPtr<phaseSystem::massTransferTable> massTransfer() const;
-
         //- Return the interfacial mass flow rate
         virtual tmp<volScalarField> iDmdt(const phasePairKey& key) const;
 
         //- Return the total interfacial mass transfer rate for phase
         virtual tmp<volScalarField> iDmdt(const phaseModel& phase) const;
 
+        //- Return the boundary mass flow rate
+        virtual tmp<volScalarField> wDmdt(const phasePairKey& key) const;
+
+        //- Return the total boundary mass transfer rate for phase
+        virtual tmp<volScalarField> wDmdt(const phaseModel& phase) const;
+
+        //- Return the interfacial mass flow rate
+        virtual tmp<volScalarField> dmdt(const phasePairKey& key) const
+        {
+            return BasePhaseSystem::dmdt(key)
+              + this->iDmdt(key) + this->wDmdt(key);
+        };
+
+        //- Return the total interfacial mass transfer rate for phase
+        virtual tmp<volScalarField> dmdt(const phaseModel& phase) const;
+
+        //- Return the momentum transfer matrices
+        virtual autoPtr<phaseSystem::momentumTransferTable>
+            momentumTransfer() const;
+
+        //- Return the heat transfer matrices
+        virtual autoPtr<phaseSystem::heatTransferTable> heatTransfer() const;
+
+        //- Return the mass transfer matrices
+        virtual autoPtr<phaseSystem::massTransferTable> massTransfer() const;
+
         //- Correct the thermodynamics
         virtual void correctThermo();
 

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/diameterModels/isothermalDiameter/isothermalDiameter.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseCompressibleTurbulenceModels/Make/options

@@ -2,8 +2,12 @@ EXE_INC = \
     -I../multiphaseSystem/lnInclude \
     -I../../phaseSystems/lnInclude \
     -I../../interfacialModels/lnInclude\
+    -I../../interfacialCompositionModels/lnInclude \
+    -I../../derivedFvPatchFields/lnInclude \
     -I$(LIB_SRC)/transportModels/compressible/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
     -I$(LIB_SRC)/transportModels/incompressible/transportModel \
     -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
     -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
@@ -16,6 +20,7 @@ LIB_LIBS = \
     -lreactingMultiphaseSystem \
     -lreactingEulerianInterfacialModels \
     -lreactingEulerianInterfacialCompositionModels \
+    -lreactingEulerianFvPatchFields \
     -lfiniteVolume \
     -lfvOptions \
     -lmeshTools

applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/Make/options

@@ -2,6 +2,8 @@ EXE_INC = \
     -I../../interfacialModels/lnInclude \
     -I../../interfacialCompositionModels/lnInclude \
     -I../../phaseSystems/lnInclude \
+    -I../multiphaseCompressibleTurbulenceModels/lnInclude \
+    -I../../derivedFvPatchFields/lnInclude \
     -IalphaContactAngle \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \

applications/solvers/multiphase/reactingEulerFoam/reactingMultiphaseEulerFoam/multiphaseSystem/multiphaseSystems.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,6 +29,7 @@ License
 #include "multiphaseSystem.H"
 #include "MomentumTransferPhaseSystem.H"
 #include "HeatTransferPhaseSystem.H"
+#include "HeatAndMassTransferPhaseSystem.H"
 #include "InterfaceCompositionPhaseChangePhaseSystem.H"
 #include "ThermalPhaseChangePhaseSystem.H"
 
@@ -52,10 +53,28 @@ namespace Foam
     );
 
     typedef
-        InterfaceCompositionPhaseChangePhaseSystem
+        HeatAndMassTransferPhaseSystem
         <
             MomentumTransferPhaseSystem<multiphaseSystem>
         >
+        heatAndMassTransferMultiphaseSystem;
+
+    addNamedToRunTimeSelectionTable
+    (
+        multiphaseSystem,
+        heatAndMassTransferMultiphaseSystem,
+        dictionary,
+        heatAndMassTransferMultiphaseSystem
+    );
+
+    typedef
+        InterfaceCompositionPhaseChangePhaseSystem
+        <
+            HeatAndMassTransferPhaseSystem
+            <
+                MomentumTransferPhaseSystem<multiphaseSystem>
+            >
+        >
         interfaceCompositionPhaseChangeMultiphaseSystem;
 
     addNamedToRunTimeSelectionTable
@@ -69,7 +88,10 @@ namespace Foam
     typedef
         ThermalPhaseChangePhaseSystem
         <
-            MomentumTransferPhaseSystem<multiphaseSystem>
+            HeatAndMassTransferPhaseSystem
+            <
+                MomentumTransferPhaseSystem<multiphaseSystem>
+            >
         >
         thermalPhaseChangeMultiphaseSystem;
 

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/Make/files

@@ -36,33 +36,5 @@ kineticTheoryModels/frictionalStressModel/JohnsonJacksonSchaeffer/JohnsonJackson
 kineticTheoryModels/derivedFvPatchFields/JohnsonJacksonParticleTheta/JohnsonJacksonParticleThetaFvPatchScalarField.C
 kineticTheoryModels/derivedFvPatchFields/JohnsonJacksonParticleSlip/JohnsonJacksonParticleSlipFvPatchVectorField.C
 
-derivedFvPatchFields/wallBoilingSubModels/partitioningModels/partitioningModel/partitioningModel.C
-derivedFvPatchFields/wallBoilingSubModels/partitioningModels/partitioningModel/newPartitioningModel.C
-derivedFvPatchFields/wallBoilingSubModels/partitioningModels/phaseFraction/phaseFraction.C
-derivedFvPatchFields/wallBoilingSubModels/partitioningModels/Lavieville/Lavieville.C
-derivedFvPatchFields/wallBoilingSubModels/partitioningModels/cosine/cosine.C
-derivedFvPatchFields/wallBoilingSubModels/partitioningModels/linear/linear.C
-
-derivedFvPatchFields/wallBoilingSubModels/nucleationSiteModels/nucleationSiteModel/nucleationSiteModel.C
-derivedFvPatchFields/wallBoilingSubModels/nucleationSiteModels/nucleationSiteModel/newNucleationSiteModel.C
-derivedFvPatchFields/wallBoilingSubModels/nucleationSiteModels/LemmertChawla/LemmertChawla.C
-
-derivedFvPatchFields/wallBoilingSubModels/departureDiameterModels/departureDiameterModel/departureDiameterModel.C
-derivedFvPatchFields/wallBoilingSubModels/departureDiameterModels/departureDiameterModel/newDepartureDiameterModel.C
-derivedFvPatchFields/wallBoilingSubModels/departureDiameterModels/TolubinskiKostanchuk/TolubinskiKostanchuk.C
-derivedFvPatchFields/wallBoilingSubModels/departureDiameterModels/KocamustafaogullariIshii/KocamustafaogullariIshii.C
-
-derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/departureFrequencyModel/departureFrequencyModel.C
-derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/departureFrequencyModel/newDepartureFrequencyModel.C
-derivedFvPatchFields/wallBoilingSubModels/departureFrequencyModels/Cole/Cole.C
-
-derivedFvPatchFields/wallBoilingSubModels/IATEsource/wallBoiling.C
-
-derivedFvPatchFields/alphatPhaseChangeJayatillekeWallFunction/alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField.C
-derivedFvPatchFields/alphatFixedDmdtWallBoilingWallFunction/alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField.C
-derivedFvPatchFields/alphatWallBoilingWallFunction/alphatWallBoilingWallFunctionFvPatchScalarField.C
-derivedFvPatchFields/copiedFixedValue/copiedFixedValueFvPatchScalarField.C
-derivedFvPatchFields/fixedMultiPhaseHeatFlux/fixedMultiPhaseHeatFluxFvPatchScalarField.C
-
 
 LIB = $(FOAM_LIBBIN)/libtwoPhaseReactingTurbulenceModels

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseCompressibleTurbulenceModels/Make/options

@@ -3,6 +3,7 @@ EXE_INC = \
     -I../../phaseSystems/lnInclude \
     -I../../interfacialModels/lnInclude\
     -I../../interfacialCompositionModels/lnInclude \
+    -I../../derivedFvPatchFields/lnInclude \
     -I$(LIB_SRC)/transportModels/compressible/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
     -I$(LIB_SRC)/transportModels/incompressible/transportModel \
@@ -17,6 +18,7 @@ LIB_LIBS = \
     -lreactingTwoPhaseSystem \
     -lreactingEulerianInterfacialModels \
     -lreactingEulerianInterfacialCompositionModels \
+    -lreactingEulerianFvPatchFields \
     -lfiniteVolume \
     -lfvOptions \
     -lmeshTools

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/Make/files

@@ -4,6 +4,7 @@ diameterModels/IATE/IATEsources/wakeEntrainmentCoalescence/wakeEntrainmentCoales
 diameterModels/IATE/IATEsources/turbulentBreakUp/turbulentBreakUp.C
 diameterModels/IATE/IATEsources/randomCoalescence/randomCoalescence.C
 diameterModels/IATE/IATEsources/phaseChange/phaseChange.C
+diameterModels/IATE/IATEsources/wallBoiling/wallBoiling.C
 
 twoPhaseSystem.C
 newTwoPhaseSystem.C

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/Make/options

@@ -2,6 +2,8 @@ EXE_INC = \
     -I../../interfacialModels/lnInclude \
     -I../../interfacialCompositionModels/lnInclude \
     -I../../phaseSystems/lnInclude \
+    -I../twoPhaseCompressibleTurbulenceModels/lnInclude \
+    -I../../derivedFvPatchFields/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/diameterModels/IATE/IATEsources/phaseChange/phaseChange.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -54,7 +54,9 @@ Foam::diameterModels::IATEsources::phaseChange::phaseChange
     const dictionary& dict
 )
 :
-    IATEsource(iate)
+    IATEsource(iate),
+    pairName_(dict.lookup("pairName")),
+    iDmdtPtr_(nullPtr)
 {}
 
 
@@ -67,24 +69,21 @@ Foam::diameterModels::IATEsources::phaseChange::R
     volScalarField& kappai
 ) const
 {
-    const ThermalPhaseChangePhaseSystem
+
-    <
+    if (!iDmdtPtr_)
-        MomentumTransferPhaseSystem
+    {
-        <
+        iDmdtPtr_ = &alphai.mesh().lookupObject<volScalarField>
-            twoPhaseSystem
+            (
-        >
+                IOobject::groupName("iDmdt",pairName_)
-    >& phaseChangeFluid = refCast
+            );
-    <
+    }
-        const ThermalPhaseChangePhaseSystem
+
-        <
+    const volScalarField& iDmdt = *iDmdtPtr_;
-            MomentumTransferPhaseSystem<twoPhaseSystem>
-        >
-    >(fluid());
 
     return -fvm::SuSp
     (
         (1.0/3.0)
-       *phaseChangeFluid.iDmdt(phase())()()
+       *iDmdt()
        /(alphai()*phase().rho()()),
         kappai
     );

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/diameterModels/IATE/IATEsources/phaseChange/phaseChange.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -54,6 +54,13 @@ class phaseChange
 :
     public IATEsource
 {
+    // Private data
+
+        //- Phase pair name
+        word pairName_;
+
+        //- Pointer to the corresponding iDmdt field
+        mutable const volScalarField* iDmdtPtr_;
 
 public:
 

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/diameterModels/IATE/IATEsources/wallBoiling/wallBoiling.C

@@ -0,0 +1,142 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "wallBoiling.H"
+#include "phaseCompressibleTurbulenceModel.H"
+#include "alphatWallBoilingWallFunctionFvPatchScalarField.H"
+#include "fvmSup.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace diameterModels
+{
+namespace IATEsources
+{
+    defineTypeNameAndDebug(wallBoiling, 0);
+    addToRunTimeSelectionTable(IATEsource, wallBoiling, dictionary);
+}
+}
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::diameterModels::IATEsources::wallBoiling::wallBoiling
+(
+    const IATE& iate,
+    const dictionary& dict
+)
+:
+    IATEsource(iate)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::tmp<Foam::fvScalarMatrix>
+Foam::diameterModels::IATEsources::wallBoiling::R
+(
+    const volScalarField& alphai,
+    volScalarField& kappai
+) const
+{
+    volScalarField::Internal R
+    (
+        IOobject
+        (
+            "wallBoiling:R",
+            phase().time().timeName(),
+            phase().mesh()
+        ),
+        phase().mesh(),
+        dimensionedScalar("R", dimless/dimTime, 0)
+    );
+
+    volScalarField::Internal Rdk
+    (
+        IOobject
+        (
+            "wallBoiling:Rdk",
+            phase().time().timeName(),
+            phase().mesh()
+        ),
+        phase().mesh(),
+        dimensionedScalar("Rdk", kappai.dimensions()/dimTime, 0)
+    );
+
+    const phaseCompressibleTurbulenceModel& turbulence =
+        phase().db().lookupObjectRef<phaseCompressibleTurbulenceModel>
+        (
+            IOobject::groupName
+            (
+                turbulenceModel::propertiesName,
+                otherPhase().name()
+            )
+        );
+
+    const tmp<volScalarField> talphat(turbulence.alphat());
+    const volScalarField::Boundary& alphatBf = talphat().boundaryField();
+
+    const scalarField& rho = phase().rho();
+
+    typedef compressible::alphatWallBoilingWallFunctionFvPatchScalarField
+        alphatWallBoilingWallFunction;
+
+    forAll(alphatBf, patchi)
+    {
+        if
+        (
+            isA<alphatWallBoilingWallFunction>(alphatBf[patchi])
+        )
+        {
+            const alphatWallBoilingWallFunction& alphatw =
+                refCast<const alphatWallBoilingWallFunction>(alphatBf[patchi]);
+
+            const scalarField& dmdt = alphatw.dmdt();
+            const scalarField& dDep = alphatw.dDeparture();
+
+            const labelList& faceCells = alphatw.patch().faceCells();
+
+            forAll(alphatw, facei)
+            {
+                if (dmdt[facei] > SMALL)
+                {
+                    const label faceCelli = faceCells[facei];
+                    R[faceCelli] =
+                        dmdt[facei]/(alphai[faceCelli]*rho[faceCelli]);
+                    Rdk[faceCelli] = R[faceCelli]*(6.0/dDep[facei]);
+                }
+            }
+        }
+    }
+
+    return Rdk - fvm::Sp(R, kappai);
+}
+
+
+// ************************************************************************* //

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/diameterModels/IATE/IATEsources/wallBoiling/wallBoiling.H

@@ -0,0 +1,98 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2016-2017 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::diameterModels::IATEsources::wallBoiling
+
+Description
+    Wall-boiling IATE source.
+
+SourceFiles
+    wallBoiling.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef wallBoiling_H
+#define wallBoiling_H
+
+#include "IATEsource.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace diameterModels
+{
+namespace IATEsources
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class wallBoiling Declaration
+\*---------------------------------------------------------------------------*/
+
+class wallBoiling
+:
+    public IATEsource
+{
+
+public:
+
+    //- Runtime type information
+    TypeName("wallBoiling");
+
+
+    // Constructors
+
+        wallBoiling
+        (
+            const IATE& iate,
+            const dictionary& dict
+        );
+
+
+    //- Destructor
+    virtual ~wallBoiling()
+    {}
+
+
+    // Member Functions
+
+        virtual tmp<fvScalarMatrix> R
+        (
+            const volScalarField& alphai,
+            volScalarField& kappai
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace IATEsources
+} // End namespace diameterModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystems.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2015-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,6 +29,7 @@ License
 #include "twoPhaseSystem.H"
 #include "MomentumTransferPhaseSystem.H"
 #include "HeatTransferPhaseSystem.H"
+#include "HeatAndMassTransferPhaseSystem.H"
 #include "InterfaceCompositionPhaseChangePhaseSystem.H"
 #include "ThermalPhaseChangePhaseSystem.H"
 
@@ -52,10 +53,28 @@ namespace Foam
     );
 
     typedef
-        InterfaceCompositionPhaseChangePhaseSystem
+        HeatAndMassTransferPhaseSystem
         <
             MomentumTransferPhaseSystem<twoPhaseSystem>
         >
+        heatAndMassTransferTwoPhaseSystem;
+
+    addNamedToRunTimeSelectionTable
+    (
+        twoPhaseSystem,
+        heatAndMassTransferTwoPhaseSystem,
+        dictionary,
+        heatAndMassTransferTwoPhaseSystem
+    );
+
+    typedef
+        InterfaceCompositionPhaseChangePhaseSystem
+        <
+            HeatAndMassTransferPhaseSystem
+            <
+                MomentumTransferPhaseSystem<twoPhaseSystem>
+            >
+        >
         interfaceCompositionPhaseChangeTwoPhaseSystem;
 
     addNamedToRunTimeSelectionTable
@@ -69,7 +88,10 @@ namespace Foam
     typedef
         ThermalPhaseChangePhaseSystem
         <
-            MomentumTransferPhaseSystem<twoPhaseSystem>
+            HeatAndMassTransferPhaseSystem
+            <
+                MomentumTransferPhaseSystem<twoPhaseSystem>
+            >
         >
         thermalPhaseChangeTwoPhaseSystem;
 

src/TurbulenceModels/compressible/Make/files

@@ -14,6 +14,5 @@ $(BCs)/convectiveHeatTransfer/convectiveHeatTransferFvPatchScalarField.C
 
 turbulentFluidThermoModels/derivedFvPatchFields/wallFunctions/alphatWallFunctions/alphatWallFunction/alphatWallFunctionFvPatchScalarField.C
 turbulentFluidThermoModels/derivedFvPatchFields/wallFunctions/alphatWallFunctions/alphatJayatillekeWallFunction/alphatJayatillekeWallFunctionFvPatchScalarField.C
-turbulentFluidThermoModels/derivedFvPatchFields/wallFunctions/alphatWallFunctions/alphatPhaseChangeWallFunction/alphatPhaseChangeWallFunctionFvPatchScalarField.C
 
 LIB = $(FOAM_LIBBIN)/libcompressibleTurbulenceModels

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/T.gas

@@ -0,0 +1,55 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "5";
+    object      T.gas;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 1 0 0 0];
+
+internalField   uniform 368;
+
+boundaryField
+{
+    inlet
+    {
+        type            fixedValue;
+        value           uniform 368;
+    }
+    outlet
+    {
+        type            inletOutlet;
+        phi             phi.gas;
+        inletValue      uniform 368;
+        value           uniform 368;
+    }
+    wall1
+    {
+        type            copiedFixedValue;
+        sourceFieldName T.liquid;
+        value           uniform 368;
+    }
+    wall2
+    {
+        type            copiedFixedValue;
+        sourceFieldName T.liquid;
+        value           uniform 368;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/T.liquid

@@ -0,0 +1,59 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "5";
+    object      T.liquid;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 1 0 0 0];
+
+internalField   uniform 368;
+
+boundaryField
+{
+    inlet
+    {
+        type            fixedValue;
+        value           uniform 368;
+    }
+    outlet
+    {
+        type            inletOutlet;
+        phi             phi.liquid;
+        inletValue      uniform 368;
+        value           uniform 368;
+    }
+    wall1
+    {
+        type            fixedMultiPhaseHeatFlux;
+        relax           0.5;
+        q               uniform 0;
+        phase           "liquid";
+        value           uniform 368;
+    }
+    wall2
+    {
+        type            fixedMultiPhaseHeatFlux;
+        relax           0.5;
+        q               uniform 0;
+        phase           "liquid";
+        value           uniform 368;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/U.gas

@@ -0,0 +1,54 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volVectorField;
+    location    "0";
+    object      U.gas;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -1 0 0 0 0];
+
+internalField   uniform (0 1 0);
+
+boundaryField
+{
+    inlet
+    {
+        type            mapped;
+        field           U.gas;
+        setAverage      1;
+        average         (0 1 0);
+        interpolationScheme cell;
+        value           uniform (0 1 0);
+    }
+    outlet
+    {
+        type            pressureInletOutletVelocity;
+        phi             phi.gas;
+        value           uniform (0 1 0);
+    }
+    wall1
+    {
+        type            slip;
+    }
+    wall2
+    {
+        type            slip;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/U.gas2

@@ -0,0 +1,54 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volVectorField;
+    location    "0";
+    object      U.gas;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -1 0 0 0 0];
+
+internalField   uniform (0 1 0);
+
+boundaryField
+{
+    inlet
+    {
+        type            mapped;
+        field           U.gas;
+        setAverage      1;
+        average         (0 1 0);
+        interpolationScheme cell;
+        value           uniform (0 1 0);
+    }
+    outlet
+    {
+        type            pressureInletOutletVelocity;
+        phi             phi.gas;
+        value           uniform (0 1 0);
+    }
+    wall1
+    {
+        type            slip;
+    }
+    wall2
+    {
+        type            slip;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/U.liquid

@@ -0,0 +1,54 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volVectorField;
+    location    "0";
+    object      U.liquid;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -1 0 0 0 0];
+
+internalField   uniform (0 1 0);
+
+boundaryField
+{
+    inlet
+    {
+        type            mapped;
+        field           U.liquid;
+        setAverage      1;
+        average         (0 1 0);
+        interpolationScheme cell;
+        value           uniform (0 1 0);
+    }
+    outlet
+    {
+        type            pressureInletOutletVelocity;
+        phi             phi.liquid;
+        value           uniform (0 1 0);
+    }
+    wall1
+    {
+        type            noSlip;
+    }
+    wall2
+    {
+        type            noSlip;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/alpha.gas

@@ -10,28 +10,28 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
-    location    "5";
+    location    "0";
-    object      water.gas;
+    object      alpha.gas;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 dimensions      [0 0 0 0 0 0 0];
 
-internalField   uniform 1;
+internalField   uniform 0;
 
 boundaryField
 {
     inlet
     {
         type            fixedValue;
-        value           uniform 1;
+        value           uniform 0;
     }
     outlet
     {
         type            inletOutlet;
         phi             phi.gas;
-        inletValue      uniform 1;
+        inletValue      uniform 0;
-        value           uniform 1;
+        value           uniform 0;
     }
     wall1
     {

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/alpha.liquid

@@ -0,0 +1,51 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      alpha.liquid;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 0 0 0 0];
+
+internalField   uniform 1;
+
+boundaryField
+{
+    inlet
+    {
+        type            calculated;
+        value           uniform 1;
+    }
+    outlet
+    {
+        type            calculated;
+        value           uniform 1;
+    }
+    wall1
+    {
+        type            calculated;
+        value           uniform 1;
+    }
+    wall2
+    {
+        type            calculated;
+        value           uniform 1;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/alphat.gas

@@ -0,0 +1,73 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      alphat.gas;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [1 -1 -1 0 0 0 0];
+
+internalField   uniform 1e-8;
+
+boundaryField
+{
+    inlet
+    {
+        type            calculated;
+        value           uniform 1e-8;
+    }
+    outlet
+    {
+        type            calculated;
+        value           uniform 1e-8;
+    }
+    wall1
+    {
+        type            compressible::alphatWallBoilingWallFunction;
+        otherPhase      liquid;
+        phaseType       vapor;
+        Prt             0.85;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        partitioningModel
+        {
+            type        phaseFraction;
+            alphaCrit   0.2;
+        }
+        value           uniform 0;
+    }
+    wall2
+    {
+        type            compressible::alphatWallBoilingWallFunction;
+        otherPhase      liquid;
+        phaseType       vapor;
+        Prt             0.85;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        partitioningModel
+        {
+            type        phaseFraction;
+            alphaCrit   0.2;
+        }
+        value           uniform 0;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/alphat.liquid

@@ -0,0 +1,101 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      alphat.liquid;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [1 -1 -1 0 0 0 0];
+
+internalField   uniform 1e-8;
+
+boundaryField
+{
+    inlet
+    {
+        type            calculated;
+        value           uniform 1e-8;
+    }
+    outlet
+    {
+        type            calculated;
+        value           uniform 1e-8;
+    }
+    wall1
+    {
+        type            compressible::alphatWallBoilingWallFunction;
+        otherPhase      gas;
+        phaseType       liquid;
+        Prt             0.85;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        relax           0.01;
+        dmdt            uniform 0;
+        partitioningModel
+        {
+            type        phaseFraction;//Lavieville;
+            alphaCrit   0.2;
+        }
+        nucleationSiteModel
+        {
+            type        LemmertChawla;
+        }
+        departureDiamModel
+        {
+            type        TolubinskiKostanchuk;
+        }
+        departureFreqModel
+        {
+            type        Cole;
+        }
+        value           uniform 0.01;
+    }
+    wall2
+    {
+        type            compressible::alphatWallBoilingWallFunction;
+        otherPhase      gas;
+        phaseType       liquid;
+        Prt             0.85;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        relax           0.01;
+        dmdt            uniform 0;
+        partitioningModel
+        {
+            type        phaseFraction;//Lavieville;
+            alphaCrit   0.2;
+        }
+        nucleationSiteModel
+        {
+            type        LemmertChawla;
+        }
+        departureDiamModel
+        {
+            type        TolubinskiKostanchuk;
+        }
+        departureFreqModel
+        {
+            type        Cole;
+        }
+        value           uniform 0.01;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/epsilon.gas

@@ -0,0 +1,62 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      epsilon.gas;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 2 -3 0 0 0 0];
+
+internalField   uniform 0.001;
+
+boundaryField
+{
+    inlet
+    {
+        type            mapped;
+        field           epsilon.gas;
+        setAverage      0;
+        average         0;
+        interpolationScheme cell;
+        value           uniform 0.01;
+    }
+    outlet
+    {
+        type            inletOutlet;
+        phi             phi.gas;
+        inletValue      uniform 0.001;
+    }
+    wall1
+    {
+        type            epsilonWallFunction;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        value           uniform 0.001;
+    }
+    wall2
+    {
+        type            epsilonWallFunction;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        value           uniform 0.001;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/epsilon.liquid

@@ -0,0 +1,63 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      epsilon.liquid;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 2 -3 0 0 0 0];
+
+internalField   uniform 0.001;
+
+boundaryField
+{
+    inlet
+    {
+        type            mapped;
+        field           epsilon.liquid;
+        setAverage      0;
+        average         0;
+        interpolationScheme cell;
+        value           uniform 0.001;
+    }
+    outlet
+    {
+        type            inletOutlet;
+        phi             phi.liquid;
+        inletValue      uniform 0.00015;
+        value           uniform 0.001;
+    }
+    wall1
+    {
+        type            epsilonWallFunction;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        value           uniform 0.001;
+    }
+    wall2
+    {
+        type            epsilonWallFunction;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        value           uniform 0.001;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/k.gas

@@ -0,0 +1,57 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      k.gas;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 2 -2 0 0 0 0];
+
+internalField   uniform 0.01;
+
+boundaryField
+{
+    inlet
+    {
+        type            mapped;
+        field           k.gas;
+        setAverage      0;
+        average         0;
+        interpolationScheme cell;
+        value           uniform 0.01;
+    }
+    outlet
+    {
+        type            inletOutlet;
+        phi             phi.gas;
+        inletValue      uniform 0.01;
+        value           uniform 0.01;
+    }
+    wall1
+    {
+        type            kqRWallFunction;
+        value           uniform 0.01;
+    }
+    wall2
+    {
+        type            kqRWallFunction;
+        value           uniform 0.01;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/k.liquid

@@ -0,0 +1,57 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      k.liquid;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 2 -2 0 0 0 0];
+
+internalField   uniform 0.01;
+
+boundaryField
+{
+    inlet
+    {
+        type            mapped;
+        field           k.liquid;
+        setAverage      0;
+        average         0;
+        interpolationScheme cell;
+        value           uniform 0.01;
+    }
+    outlet
+    {
+        type            inletOutlet;
+        phi             phi.liquid;
+        inletValue      uniform 3.75e-05;
+        value           uniform 0.01;
+    }
+    wall1
+    {
+        type            kqRWallFunction;
+        value           uniform 0.01;
+    }
+    wall2
+    {
+        type            kqRWallFunction;
+        value           uniform 0.01;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/reactingMultiphaseEulerFoam/RAS/wallBoiling1D_2phase/0.org/nut.gas

@@ -0,0 +1,57 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      nut.gas;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 2 -1 0 0 0 0];
+
+internalField   uniform 0.0001;
+
+boundaryField
+{
+    inlet
+    {
+        type            calculated;
+        value           uniform 0.0001;
+    }
+    outlet
+    {
+        type            calculated;
+        value           uniform 0.0001;
+    }
+    wall1
+    {
+        type            nutkWallFunction;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        value           uniform 0;
+    }
+    wall2
+    {
+        type            nutkWallFunction;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        value           uniform 0;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //