ENH: Updated regionModels from internal development line

src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -52,8 +52,10 @@ void noPyrolysis::constructThermoChemistry()
         basicSolidChemistryModel::New(regionMesh()).ptr()
     );
 
-    solidThermo_.reset(&solidChemistry_->solidThermo());
+    radiation_.reset(radiation::radiationModel::New
-    radiation_.reset(radiation::radiationModel::New(solidThermo_->T()).ptr());
+    (
+        solidChemistry_->thermo().T()
+    ).ptr());
 }
 
 
@@ -96,7 +98,6 @@ noPyrolysis::noPyrolysis
 :
     pyrolysisModel(mesh, regionType),
     solidChemistry_(NULL),
-    solidThermo_(NULL),
     radiation_(NULL)
 {
     if (active())
@@ -116,7 +117,6 @@ noPyrolysis::noPyrolysis
 :
     pyrolysisModel(mesh, regionType),
     solidChemistry_(NULL),
-    solidThermo_(NULL),
     radiation_(NULL)
 {
     if (active())
@@ -148,19 +148,19 @@ void noPyrolysis::evolveRegion()
 
 const volScalarField& noPyrolysis::rho() const
 {
-    return solidThermo_->rho();
+    return solidChemistry_->thermo().rho();
 }
 
 
 const volScalarField& noPyrolysis::T() const
 {
-    return solidThermo_->T();
+    return solidChemistry_->thermo().T();
 }
 
 
 const tmp<volScalarField> noPyrolysis::Cp() const
 {
-    return solidThermo_->Cp();
+    return solidChemistry_->thermo().Cp();
 }
 
 
@@ -172,7 +172,7 @@ tmp<volScalarField> noPyrolysis::kappaRad() const
 
 tmp<volScalarField> noPyrolysis::kappa() const
 {
-     return solidThermo_->kappa();
+    return solidChemistry_->thermo().kappa();
 }
 
 

src/regionModels/pyrolysisModels/noPyrolysis/noPyrolysis.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -81,12 +81,9 @@ protected:
         //- Reset solidChemistryModel and solidThermo pointers
         void constructThermoChemistry();
 
-        //- Reference to the solid chemistry model
+        //- Pointer to the solid chemistry model
         autoPtr<basicSolidChemistryModel> solidChemistry_;
 
-        //- Reference to solid thermo
-        autoPtr<solidReactionThermo> solidThermo_;
-
         //- Pointer to radiation model
         autoPtr<radiation::radiationModel> radiation_;
 

src/regionModels/pyrolysisModels/pyrolysisModel/pyrolysisModelCollection.C

@@ -53,7 +53,6 @@ namespace pyrolysisModels
 pyrolysisModelCollection::pyrolysisModelCollection(const fvMesh& mesh)
 :
     PtrList<pyrolysisModel>()
-
 {
     IOdictionary pyrolysisZonesDict
     (
@@ -162,12 +161,9 @@ scalar pyrolysisModelCollection::maxDiff() const
     scalar maxDiff = 0.0;
     forAll(*this, i)
     {
-        if (maxDiff < this->operator[](i).maxDiff())
+        maxDiff = max(maxDiff, this->operator[](i).maxDiff());
-        {
-            maxDiff = this->operator[](i).maxDiff();
     }
 
-    }
     return maxDiff;
 }
 
@@ -175,16 +171,9 @@ scalar pyrolysisModelCollection::maxDiff() const
 scalar pyrolysisModelCollection::solidRegionDiffNo() const
 {
     scalar totalDiNum = GREAT;
-
     forAll(*this, i)
     {
-        if
+        totalDiNum = min(totalDiNum, this->operator[](i).solidRegionDiffNo());
-        (
-            totalDiNum > this->operator[](i).solidRegionDiffNo()
-        )
-        {
-            totalDiNum = this->operator[](i).solidRegionDiffNo();
-        }
     }
 
     return totalDiNum;

src/regionModels/surfaceFilmModels/derivedFvPatchFields/inclinedFilmNusseltHeight/inclinedFilmNusseltHeightFvPatchScalarField.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -38,6 +38,7 @@ inclinedFilmNusseltHeightFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(p, iF),
+    filmRegionName_("surfaceFilmProperties"),
     GammaMean_(),
     a_(),
     omega_()
@@ -54,6 +55,7 @@ inclinedFilmNusseltHeightFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(ptf, p, iF, mapper),
+    filmRegionName_(ptf.filmRegionName_),
     GammaMean_(ptf.GammaMean_().clone().ptr()),
     a_(ptf.a_().clone().ptr()),
     omega_(ptf.omega_().clone().ptr())
@@ -69,6 +71,10 @@ inclinedFilmNusseltHeightFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(p, iF),
+    filmRegionName_
+    (
+        dict.lookupOrDefault<word>("filmRegion", "surfaceFilmProperties")
+    ),
     GammaMean_(DataEntry<scalar>::New("GammaMean", dict)),
     a_(DataEntry<scalar>::New("a", dict)),
     omega_(DataEntry<scalar>::New("omega", dict))
@@ -84,6 +90,7 @@ inclinedFilmNusseltHeightFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(wmfrhpsf),
+    filmRegionName_(wmfrhpsf.filmRegionName_),
     GammaMean_(wmfrhpsf.GammaMean_().clone().ptr()),
     a_(wmfrhpsf.a_().clone().ptr()),
     omega_(wmfrhpsf.omega_().clone().ptr())
@@ -98,6 +105,7 @@ inclinedFilmNusseltHeightFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(wmfrhpsf, iF),
+    filmRegionName_(wmfrhpsf.filmRegionName_),
     GammaMean_(wmfrhpsf.GammaMean_().clone().ptr()),
     a_(wmfrhpsf.a_().clone().ptr()),
     omega_(wmfrhpsf.omega_().clone().ptr())
@@ -118,10 +126,7 @@ void Foam::inclinedFilmNusseltHeightFvPatchScalarField::updateCoeffs()
     // retrieve the film region from the database
 
     const regionModels::regionModel& region =
-        db().time().lookupObject<regionModels::regionModel>
+        db().time().lookupObject<regionModels::regionModel>(filmRegionName_);
-        (
-            "surfaceFilmProperties"
-        );
 
     const regionModels::surfaceFilmModels::kinematicSingleLayer& film =
         dynamic_cast
@@ -134,8 +139,7 @@ void Foam::inclinedFilmNusseltHeightFvPatchScalarField::updateCoeffs()
     // note: normal pointing into the domain
     const vectorField n(-patch().nf());
 
-    // TODO: currently re-evaluating the entire gTan field to return this patch
+    const scalarField gTan(film.gTan(patchI) & n);
-    const scalarField gTan(film.gTan()().boundaryField()[patchI] & n);
 
     if (patch().size() && (max(mag(gTan)) < SMALL))
     {
@@ -195,6 +199,13 @@ void Foam::inclinedFilmNusseltHeightFvPatchScalarField::write
 ) const
 {
     fixedValueFvPatchScalarField::write(os);
+    writeEntryIfDifferent<word>
+    (
+        os,
+        "filmRegion",
+        "surfaceFilmProperties",
+        filmRegionName_
+    );
     GammaMean_->writeData(os);
     a_->writeData(os);
     omega_->writeData(os);

src/regionModels/surfaceFilmModels/derivedFvPatchFields/inclinedFilmNusseltHeight/inclinedFilmNusseltHeightFvPatchScalarField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -56,6 +56,9 @@ class inclinedFilmNusseltHeightFvPatchScalarField
 {
     // Private data
 
+        //- Name of film region
+        word filmRegionName_;
+
         //- Mean mass flow rate per unit length [kg/s/m]
         autoPtr<DataEntry<scalar> > GammaMean_;
 

src/regionModels/surfaceFilmModels/derivedFvPatchFields/inclinedFilmNusseltInletVelocity/inclinedFilmNusseltInletVelocityFvPatchVectorField.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -38,6 +38,7 @@ inclinedFilmNusseltInletVelocityFvPatchVectorField
 )
 :
     fixedValueFvPatchVectorField(p, iF),
+    filmRegionName_("surfaceFilmProperties"),
     GammaMean_(),
     a_(),
     omega_()
@@ -54,6 +55,7 @@ inclinedFilmNusseltInletVelocityFvPatchVectorField
 )
 :
     fixedValueFvPatchVectorField(ptf, p, iF, mapper),
+    filmRegionName_(ptf.filmRegionName_),
     GammaMean_(ptf.GammaMean_().clone().ptr()),
     a_(ptf.a_().clone().ptr()),
     omega_(ptf.omega_().clone().ptr())
@@ -69,6 +71,10 @@ inclinedFilmNusseltInletVelocityFvPatchVectorField
 )
 :
     fixedValueFvPatchVectorField(p, iF),
+    filmRegionName_
+    (
+        dict.lookupOrDefault<word>("filmRegion", "surfaceFilmProperties")
+    ),
     GammaMean_(DataEntry<scalar>::New("GammaMean", dict)),
     a_(DataEntry<scalar>::New("a", dict)),
     omega_(DataEntry<scalar>::New("omega", dict))
@@ -84,6 +90,7 @@ inclinedFilmNusseltInletVelocityFvPatchVectorField
 )
 :
     fixedValueFvPatchVectorField(fmfrpvf),
+    filmRegionName_(fmfrpvf.filmRegionName_),
     GammaMean_(fmfrpvf.GammaMean_().clone().ptr()),
     a_(fmfrpvf.a_().clone().ptr()),
     omega_(fmfrpvf.omega_().clone().ptr())
@@ -98,6 +105,7 @@ inclinedFilmNusseltInletVelocityFvPatchVectorField
 )
 :
     fixedValueFvPatchVectorField(fmfrpvf, iF),
+    filmRegionName_(fmfrpvf.filmRegionName_),
     GammaMean_(fmfrpvf.GammaMean_().clone().ptr()),
     a_(fmfrpvf.a_().clone().ptr()),
     omega_(fmfrpvf.omega_().clone().ptr())
@@ -118,10 +126,7 @@ void Foam::inclinedFilmNusseltInletVelocityFvPatchVectorField::updateCoeffs()
     // retrieve the film region from the database
 
     const regionModels::regionModel& region =
-        db().time().lookupObject<regionModels::regionModel>
+        db().time().lookupObject<regionModels::regionModel>(filmRegionName_);
-        (
-            "surfaceFilmProperties"
-        );
 
     const regionModels::surfaceFilmModels::kinematicSingleLayer& film =
         dynamic_cast
@@ -133,8 +138,7 @@ void Foam::inclinedFilmNusseltInletVelocityFvPatchVectorField::updateCoeffs()
     // note: normal pointing into the domain
     const vectorField n(-patch().nf());
 
-    // TODO: currently re-evaluating the entire gTan field to return this patch
+    const scalarField gTan(film.gTan(patchI) & n);
-    const scalarField gTan(film.gTan()().boundaryField()[patchI] & n);
 
     if (patch().size() && (max(mag(gTan)) < SMALL))
     {
@@ -191,6 +195,13 @@ void Foam::inclinedFilmNusseltInletVelocityFvPatchVectorField::write
 ) const
 {
     fvPatchVectorField::write(os);
+    writeEntryIfDifferent<word>
+    (
+        os,
+        "filmRegion",
+        "surfaceFilmProperties",
+        filmRegionName_
+    );
     GammaMean_->writeData(os);
     a_->writeData(os);
     omega_->writeData(os);

src/regionModels/surfaceFilmModels/derivedFvPatchFields/inclinedFilmNusseltInletVelocity/inclinedFilmNusseltInletVelocityFvPatchVectorField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -56,6 +56,9 @@ class inclinedFilmNusseltInletVelocityFvPatchVectorField
 {
     // Private data
 
+        //- Name of film region
+        word filmRegionName_;
+
         //- Mean mass flow rate per unit length [kg/s/m]
         autoPtr<DataEntry<scalar> > GammaMean_;
 

src/regionModels/surfaceFilmModels/derivedFvPatchFields/wallFunctions/alphatFilmWallFunction/alphatFilmWallFunctionFvPatchScalarField.C

@@ -51,6 +51,7 @@ alphatFilmWallFunctionFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(p, iF),
+    filmRegionName_("surfaceFilmProperties"),
     B_(5.5),
     yPlusCrit_(11.05),
     Cmu_(0.09),
@@ -69,6 +70,7 @@ alphatFilmWallFunctionFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(ptf, p, iF, mapper),
+    filmRegionName_(ptf.filmRegionName_),
     B_(ptf.B_),
     yPlusCrit_(ptf.yPlusCrit_),
     Cmu_(ptf.Cmu_),
@@ -86,6 +88,10 @@ alphatFilmWallFunctionFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(p, iF, dict),
+    filmRegionName_
+    (
+        dict.lookupOrDefault<word>("filmRegion", "surfaceFilmProperties")
+    ),
     B_(dict.lookupOrDefault("B", 5.5)),
     yPlusCrit_(dict.lookupOrDefault("yPlusCrit", 11.05)),
     Cmu_(dict.lookupOrDefault("Cmu", 0.09)),
@@ -101,6 +107,7 @@ alphatFilmWallFunctionFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(fwfpsf),
+    filmRegionName_(fwfpsf.filmRegionName_),
     B_(fwfpsf.B_),
     yPlusCrit_(fwfpsf.yPlusCrit_),
     Cmu_(fwfpsf.Cmu_),
@@ -117,6 +124,7 @@ alphatFilmWallFunctionFvPatchScalarField
 )
 :
     fixedValueFvPatchScalarField(fwfpsf, iF),
+    filmRegionName_(fwfpsf.filmRegionName_),
     B_(fwfpsf.B_),
     yPlusCrit_(fwfpsf.yPlusCrit_),
     Cmu_(fwfpsf.Cmu_),
@@ -141,8 +149,7 @@ void alphatFilmWallFunctionFvPatchScalarField::updateCoeffs()
     int oldTag = UPstream::msgType();
     UPstream::msgType() = oldTag+1;
 
-    bool foundFilm =
+    bool foundFilm = db().time().foundObject<modelType>(filmRegionName_);
-        db().time().foundObject<modelType>("surfaceFilmProperties");
 
     if (!foundFilm)
     {
@@ -154,7 +161,7 @@ void alphatFilmWallFunctionFvPatchScalarField::updateCoeffs()
 
     // Retrieve phase change mass from surface film model
     const modelType& filmModel =
-        db().time().lookupObject<modelType>("surfaceFilmProperties");
+        db().time().lookupObject<modelType>(filmRegionName_);
 
     const label filmPatchI = filmModel.regionPatchID(patchi);
 
@@ -230,6 +237,13 @@ void alphatFilmWallFunctionFvPatchScalarField::updateCoeffs()
 void alphatFilmWallFunctionFvPatchScalarField::write(Ostream& os) const
 {
     fvPatchField<scalar>::write(os);
+    writeEntryIfDifferent<word>
+    (
+        os,
+        "filmRegion",
+        "surfaceFilmProperties",
+        filmRegionName_
+    );
     os.writeKeyword("B") << B_ << token::END_STATEMENT << nl;
     os.writeKeyword("yPlusCrit") << yPlusCrit_ << token::END_STATEMENT << nl;
     os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;

src/regionModels/surfaceFilmModels/derivedFvPatchFields/wallFunctions/alphatFilmWallFunction/alphatFilmWallFunctionFvPatchScalarField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -92,6 +92,9 @@ protected:
 
     // Protected data
 
+        //- Name of film region
+        word filmRegionName_;
+
         //- B Coefficient (default = 5.5)
         scalar B_;
 

src/regionModels/surfaceFilmModels/derivedFvPatchFields/wallFunctions/nutkFilmWallFunction/nutkFilmWallFunctionFvPatchScalarField.C

@@ -53,8 +53,7 @@ tmp<scalarField> nutkFilmWallFunctionFvPatchScalarField::calcUTau
 
     typedef regionModels::surfaceFilmModels::surfaceFilmModel modelType;
 
-    bool foundFilm =
+    bool foundFilm = db().time().foundObject<modelType>(filmRegionName_);
-        db().time().foundObject<modelType>("surfaceFilmProperties");
 
     if (!foundFilm)
     {
@@ -66,7 +65,7 @@ tmp<scalarField> nutkFilmWallFunctionFvPatchScalarField::calcUTau
 
     // Retrieve phase change mass from surface film model
     const modelType& filmModel =
-        db().time().lookupObject<modelType>("surfaceFilmProperties");
+        db().time().lookupObject<modelType>(filmRegionName_);
 
     const label filmPatchI = filmModel.regionPatchID(patchi);
 
@@ -158,6 +157,7 @@ nutkFilmWallFunctionFvPatchScalarField::nutkFilmWallFunctionFvPatchScalarField
 )
 :
     nutkWallFunctionFvPatchScalarField(p, iF),
+    filmRegionName_("surfaceFilmProperties"),
     B_(5.5),
     yPlusCrit_(11.05)
 {}
@@ -172,6 +172,7 @@ nutkFilmWallFunctionFvPatchScalarField::nutkFilmWallFunctionFvPatchScalarField
 )
 :
     nutkWallFunctionFvPatchScalarField(ptf, p, iF, mapper),
+    filmRegionName_(ptf.filmRegionName_),
     B_(5.5),
     yPlusCrit_(11.05)
 {}
@@ -185,6 +186,10 @@ nutkFilmWallFunctionFvPatchScalarField::nutkFilmWallFunctionFvPatchScalarField
 )
 :
     nutkWallFunctionFvPatchScalarField(p, iF, dict),
+    filmRegionName_
+    (
+        dict.lookupOrDefault<word>("filmRegion", "surfaceFilmProperties")
+    ),
     B_(dict.lookupOrDefault("B", 5.5)),
     yPlusCrit_(dict.lookupOrDefault("yPlusCrit", 11.05))
 {}
@@ -196,6 +201,7 @@ nutkFilmWallFunctionFvPatchScalarField::nutkFilmWallFunctionFvPatchScalarField
 )
 :
     nutkWallFunctionFvPatchScalarField(wfpsf),
+    filmRegionName_(wfpsf.filmRegionName_),
     B_(wfpsf.B_),
     yPlusCrit_(wfpsf.yPlusCrit_)
 {}
@@ -208,6 +214,7 @@ nutkFilmWallFunctionFvPatchScalarField::nutkFilmWallFunctionFvPatchScalarField
 )
 :
     nutkWallFunctionFvPatchScalarField(wfpsf, iF),
+    filmRegionName_(wfpsf.filmRegionName_),
     B_(wfpsf.B_),
     yPlusCrit_(wfpsf.yPlusCrit_)
 {}
@@ -240,6 +247,13 @@ tmp<scalarField> nutkFilmWallFunctionFvPatchScalarField::yPlus() const
 void nutkFilmWallFunctionFvPatchScalarField::write(Ostream& os) const
 {
     fvPatchField<scalar>::write(os);
+    writeEntryIfDifferent<word>
+    (
+        os,
+        "filmRegion",
+        "surfaceFilmProperties",
+        filmRegionName_
+    );
     writeLocalEntries(os);
     os.writeKeyword("B") << B_ << token::END_STATEMENT << nl;
     os.writeKeyword("yPlusCrit") << yPlusCrit_ << token::END_STATEMENT << nl;

src/regionModels/surfaceFilmModels/derivedFvPatchFields/wallFunctions/nutkFilmWallFunction/nutkFilmWallFunctionFvPatchScalarField.H

@@ -51,8 +51,8 @@ SourceFiles
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef compressibleNutSpalartAllmarasWallFunctionFvPatchScalarField_H
+#ifndef nutkFilmWallFunctionFvPatchScalarField_H
-#define compressibleNutSpalartAllmarasWallFunctionFvPatchScalarField_H
+#define nutkFilmWallFunctionFvPatchScalarField_H
 
 #include "nutkWallFunctionFvPatchScalarField.H"
 
@@ -77,6 +77,9 @@ protected:
 
     // Protected data
 
+        //- Name of film region
+        word filmRegionName_;
+
         //- B Coefficient (default = 5.5)
         scalar B_;
 

src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.C

@@ -140,6 +140,13 @@ void kinematicSingleLayer::transferPrimaryRegionSourceFields()
     // update addedMassTotal counter
     if (time().outputTime())
     {
+        if (debug)
+        {
+            rhoSp_.write();
+            USp_.write();
+            pSp_.write();
+        }
+
         scalar addedMassTotal = 0.0;
         outputProperties().readIfPresent("addedMassTotal", addedMassTotal);
         addedMassTotal += returnReduce(addedMassTotal_, sumOp<scalar>());
@@ -422,6 +429,9 @@ void kinematicSingleLayer::solveThickness
 
     U_.correctBoundaryConditions();
 
+    // Update film wall and surface velocities
+    updateSurfaceVelocities();
+
     // Continuity check
     continuityCheck();
 }
@@ -865,6 +875,10 @@ void kinematicSingleLayer::preEvolveRegion()
 
     transferPrimaryRegionSourceFields();
 
+    updateSurfaceVelocities();
+
+    correctAlpha();
+
     // Reset transfer fields
 //    availableMass_ = mass();
     availableMass_ = netMass();
@@ -880,12 +894,6 @@ void kinematicSingleLayer::evolveRegion()
         Info<< "kinematicSingleLayer::evolveRegion()" << endl;
     }
 
-    // Update film coverage indicator
-    correctAlpha();
-
-    // Update film wall and surface velocities
-    updateSurfaceVelocities();
-
     // Update sub-models to provide updated source contributions
     updateSubmodels();
 
@@ -913,6 +921,15 @@ void kinematicSingleLayer::evolveRegion()
 
     // Update deltaRho_ with new delta_
     deltaRho_ == delta_*rho_;
+}
+
+
+void kinematicSingleLayer::postEvolveRegion()
+{
+    if (debug)
+    {
+        Info<< "kinematicSingleLayer::postEvolveRegion()" << endl;
+    }
 
     // Reset source terms for next time integration
     resetPrimaryRegionSourceTerms();
@@ -933,7 +950,7 @@ scalar kinematicSingleLayer::CourantNumber() const
 
         forAll(delta_, i)
         {
-            if (delta_[i] > deltaCoLimit_)
+            if ((delta_[i] > deltaCoLimit_) && (alpha_[i] > 0.5))
             {
                 CoNum = max(CoNum, sumPhi[i]/(delta_[i]*magSf()[i]));
             }

src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -495,6 +495,9 @@ public:
             //- Return the gravity tangential component contributions
             inline tmp<volVectorField> gTan() const;
 
+            //- Return the gravity tangential component contributions for patchI
+            inline tmp<vectorField> gTan(const label patchI) const;
+
 
         // Evolution
 
@@ -504,6 +507,9 @@ public:
             //- Evolve the film equations
             virtual void evolveRegion();
 
+            //- Post-evolve film hook
+            virtual void postEvolveRegion();
+
 
         // Source fields
 

src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayerI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -190,7 +190,7 @@ inline const filmTurbulenceModel& kinematicSingleLayer::turbulence() const
 
 inline tmp<volScalarField> kinematicSingleLayer::mass() const
 {
-    return rho_*delta_*magSf();
+    return deltaRho_*magSf();
 }
 
 
@@ -199,7 +199,7 @@ inline tmp<volScalarField> kinematicSingleLayer::netMass() const
     return
         fvc::surfaceSum(pos(phi_)*phi_/(fvc::interpolate(delta_) + deltaSmall_))
        *time().deltaT()
-      + rho_*delta_*magSf();
+      + deltaRho_*magSf();
 }
 
 
@@ -280,6 +280,18 @@ inline tmp<volVectorField> kinematicSingleLayer::gTan() const
     return tgTan;
 }
 
+inline tmp<vectorField> kinematicSingleLayer::gTan
+(
+    const label patchI
+) const
+{
+    const vectorField& nH = nHat().boundaryField()[patchI];
+    const vector& g = g_.value();
+    tmp<vectorField> tgTan(new vectorField(g - nH*(g & nH)));
+
+    return tgTan;
+}
+
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 

src/regionModels/surfaceFilmModels/submodels/kinematic/filmThermoModel/liquidFilmThermo/liquidFilmThermo.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -261,10 +261,7 @@ tmp<volScalarField> liquidFilmThermo::rho() const
 
     if (useReferenceValues_)
     {
-        forAll(rho, cellI)
+        rho = this->rho(pRef_, TRef_);
-        {
-            rho[cellI] = this->rho(pRef_, TRef_);
-        }
     }
     else
     {
@@ -309,10 +306,7 @@ tmp<volScalarField> liquidFilmThermo::mu() const
 
     if (useReferenceValues_)
     {
-        forAll(mu, cellI)
+        mu = this->mu(pRef_, TRef_);
-        {
-            mu[cellI] = this->mu(pRef_, TRef_);
-        }
     }
     else
     {
@@ -357,10 +351,7 @@ tmp<volScalarField> liquidFilmThermo::sigma() const
 
     if (useReferenceValues_)
     {
-        forAll(sigma, cellI)
+        sigma = this->sigma(pRef_, TRef_);
-        {
-            sigma[cellI] = this->sigma(pRef_, TRef_);
-        }
     }
     else
     {
@@ -405,10 +396,7 @@ tmp<volScalarField> liquidFilmThermo::Cp() const
 
     if (useReferenceValues_)
     {
-        forAll(Cp, cellI)
+        Cp = this->Cp(pRef_, TRef_);
-        {
-            Cp[cellI] = this->Cp(pRef_, TRef_);
-        }
     }
     else
     {
@@ -453,10 +441,7 @@ tmp<volScalarField> liquidFilmThermo::kappa() const
 
     if (useReferenceValues_)
     {
-        forAll(kappa, cellI)
+        kappa = this->kappa(pRef_, TRef_);
-        {
-            kappa[cellI] = this->kappa(pRef_, TRef_);
-        }
     }
     else
     {

src/regionModels/surfaceFilmModels/submodels/kinematic/force/thermocapillaryForce/thermocapillaryForce.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -63,12 +63,13 @@ thermocapillaryForce::~thermocapillaryForce()
 
 tmp<fvVectorMatrix> thermocapillaryForce::correct(volVectorField& U)
 {
+    const volScalarField& alpha = owner_.alpha();
     const volScalarField& sigma = owner_.sigma();
 
     tmp<fvVectorMatrix>
         tfvm(new fvVectorMatrix(U, dimForce/dimArea*dimVolume));
 
-    tfvm() += fvc::grad(sigma);
+    tfvm() += alpha*fvc::grad(sigma);
 
     return tfvm;
 }

src/regionModels/surfaceFilmModels/thermoSingleLayer/thermoSingleLayer.C

@@ -27,6 +27,7 @@ License
 #include "fvcDiv.H"
 #include "fvcLaplacian.H"
 #include "fvm.H"
+#include "fvcDdt.H"
 #include "addToRunTimeSelectionTable.H"
 #include "zeroGradientFvPatchFields.H"
 #include "mappedFieldFvPatchField.H"
@@ -188,6 +189,14 @@ void thermoSingleLayer::transferPrimaryRegionSourceFields()
 
     // Apply enthalpy source as difference between incoming and actual states
     hsSp_ -= rhoSp_*hs_;
+
+    if (time().outputTime())
+    {
+        if (debug)
+        {
+            hsSp_.write();
+        }
+    }
 }
 
 
@@ -231,10 +240,12 @@ void thermoSingleLayer::updateSubmodels()
     htcs_->correct();
     htcw_->correct();
 
+    scalarField availableMass(alpha_*availableMass_);
+
     phaseChange_->correct
     (
         time_.deltaTValue(),
-        availableMass_,
+        availableMass,
         primaryMassPCTrans_,
         primaryEnergyPCTrans_
     );
@@ -256,6 +267,12 @@ void thermoSingleLayer::updateSubmodels()
 
 tmp<fvScalarMatrix> thermoSingleLayer::q(volScalarField& hs) const
 {
+//    Only apply heat transfer where the film is present
+//    - leads to temperature unboundedness?
+//    volScalarField boundedAlpha(max(alpha_, ROOTVSMALL));
+//    volScalarField htcst(htcs_->h()*boundedAlpha);
+//    volScalarField htcwt(htcw_->h()*boundedAlpha);
+
     return
     (
       - fvm::Sp(htcs_->h()/Cp_, hs)
@@ -274,7 +291,13 @@ void thermoSingleLayer::solveEnergy()
         Info<< "thermoSingleLayer::solveEnergy()" << endl;
     }
 
-    updateSurfaceTemperatures();
+
+    dimensionedScalar residualDeltaRho
+    (
+        "residualDeltaRho",
+        deltaRho_.dimensions(),
+        1e-10
+    );
 
     solve
     (
@@ -282,16 +305,21 @@ void thermoSingleLayer::solveEnergy()
       + fvm::div(phi_, hs_)
      ==
       - hsSp_
-      + q(hs_)
+      + alpha_*q(hs_)
-      + radiation_->Shs()
+      + alpha_*radiation_->Shs()
 //      - fvm::SuSp(rhoSp_, hs_)
       - rhoSp_*hs_
+      + fvc::ddt(residualDeltaRho + deltaRho_, hs_)
+      - fvm::ddt(residualDeltaRho + deltaRho_, hs_)
     );
 
     correctThermoFields();
 
     // evaluate viscosity from user-model
     viscosity_->correct(pPrimary_, T_);
+
+    // Update film wall and surface temperatures
+    updateSurfaceTemperatures();
 }
 
 
@@ -611,10 +639,10 @@ void thermoSingleLayer::preEvolveRegion()
         Info<< "thermoSingleLayer::preEvolveRegion()" << endl;
     }
 
-//    correctHsForMappedT();
-
     kinematicSingleLayer::preEvolveRegion();
 
+    updateSurfaceTemperatures();
+
     // Update phase change
     primaryMassPCTrans_ == dimensionedScalar("zero", dimMass, 0.0);
     primaryEnergyPCTrans_ == dimensionedScalar("zero", dimEnergy, 0.0);
@@ -628,15 +656,6 @@ void thermoSingleLayer::evolveRegion()
         Info<< "thermoSingleLayer::evolveRegion()" << endl;
     }
 
-    // Update film coverage indicator
-    correctAlpha();
-
-    // Update film wall and surface velocities
-    updateSurfaceVelocities();
-
-    // Update film wall and surface temperatures
-    updateSurfaceTemperatures();
-
     // Update sub-models to provide updated source contributions
     updateSubmodels();
 
@@ -670,9 +689,6 @@ void thermoSingleLayer::evolveRegion()
 
     // Update temperature using latest hs_
     T_ == T(hs_);
-
-    // Reset source terms for next time integration
-    resetPrimaryRegionSourceTerms();
 }
 
 
@@ -739,7 +755,7 @@ tmp<DimensionedField<scalar, volMesh> > thermoSingleLayer::Srho() const
         (
             IOobject
             (
-                "thermoSingleLayer::Srho",
+                typeName + ":Srho",
                 time().timeName(),
                 primaryMesh(),
                 IOobject::NO_READ,
@@ -791,7 +807,7 @@ tmp<DimensionedField<scalar, volMesh> > thermoSingleLayer::Srho
         (
             IOobject
             (
-                "thermoSingleLayer::Srho(" + Foam::name(i) + ")",
+                typeName + ":Srho(" + Foam::name(i) + ")",
                 time_.timeName(),
                 primaryMesh(),
                 IOobject::NO_READ,
@@ -841,7 +857,7 @@ tmp<DimensionedField<scalar, volMesh> > thermoSingleLayer::Sh() const
         (
             IOobject
             (
-                "thermoSingleLayer::Sh",
+                typeName + ":Sh",
                 time().timeName(),
                 primaryMesh(),
                 IOobject::NO_READ,

tutorials/lagrangian/reactingParcelFilmFoam/cylinder/system/wallFilmRegion/fvSchemes

@@ -17,11 +17,7 @@ FoamFile
 
 ddtSchemes
 {
-    default none;
+    default         Euler;
-    ddt(deltaf*rhof) Euler;
-    ddt(rhof,deltaf) Euler;
-    ddt(deltaf*rhof,Uf) Euler;
-    ddt(deltaf*rhof,hf) Euler;
 }
 
 divSchemes

tutorials/lagrangian/reactingParcelFilmFoam/hotBoxes/system/wallFilmRegion.org/fvSchemes

@@ -17,11 +17,7 @@ FoamFile
 
 ddtSchemes
 {
-    default none;
+    default         Euler;
-    ddt(deltaf*rhof) Euler;
-    ddt(rhof,deltaf) Euler;
-    ddt(deltaf*rhof,Uf) Euler;
-    ddt(deltaf*rhof,hf) Euler;
 }
 
 divSchemes

tutorials/lagrangian/reactingParcelFilmFoam/rivuletPanel/system/wallFilmRegion/fvSchemes

@@ -17,11 +17,7 @@ FoamFile
 
 ddtSchemes
 {
-    default none;
+    default         Euler;
-    ddt(deltaf*rhof) Euler;
-    ddt(rhof,deltaf) Euler;
-    ddt(deltaf*rhof,Uf) Euler;
-    ddt(deltaf*rhof,hf) Euler;
 }
 
 divSchemes

tutorials/lagrangian/reactingParcelFilmFoam/splashPanel/system/wallFilmRegion/fvSchemes

@@ -17,11 +17,7 @@ FoamFile
 
 ddtSchemes
 {
-    default none;
+    default         Euler;
-    ddt(deltaf*rhof) Euler;
-    ddt(rhof,deltaf) Euler;
-    ddt(deltaf*rhof,Uf) Euler;
-    ddt(deltaf*rhof,hf) Euler;
 }
 
 divSchemes