Commit 5bab2879 authored by Henry Weller's avatar Henry Weller Committed by Andrew Heather
Browse files

compressibleInterFilmFoam: Experimental VoF solver supporting VoF<->film transfer

parent cc9ffdff
......@@ -5,5 +5,6 @@ wclean libso twoPhaseMixtureThermo
wclean libso surfaceTensionModels
wclean
wclean compressibleInterDyMFoam
wclean compressibleInterFilmFoam
#------------------------------------------------------------------------------
......@@ -9,5 +9,6 @@ wmake $targetType surfaceTensionModels
wmake $targetType
wmake $targetType compressibleInterDyMFoam
wmake $targetType compressibleInterFilmFoam
#------------------------------------------------------------------------------
VoFPatchTransfer/VoFPatchTransfer.C
compressibleInterFilmFoam.C
EXE = $(FOAM_APPBIN)/compressibleInterFilmFoam
EXE_INC = \
-I. \
-I.. \
-I../../VoF \
-I../twoPhaseMixtureThermo \
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \
-I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/thermophysicalFunctions/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
EXE_LIBS = \
-ltwoPhaseMixtureThermo \
-ltwoPhaseSurfaceTension \
-lcompressibleTransportModels \
-lfluidThermophysicalModels \
-lspecie \
-ltwoPhaseMixture \
-ltwoPhaseProperties \
-linterfaceProperties \
-lturbulenceModels \
-lcompressibleTurbulenceModels \
-lSLGThermo \
-lsurfaceFilmModels \
-lsurfaceFilmDerivedFvPatchFields \
-ldynamicMesh \
-lmeshTools \
-ldynamicFvMesh \
-lfiniteVolume \
-lfvOptions
{
fvScalarMatrix TEqn
(
fvm::ddt(rho, T) + fvm::div(rhoPhi, T)
- fvm::Sp(contErr, T)
- fvm::laplacian(mixture.alphaEff(turbulence->mut()), T)
+ (
(
fvc::div(fvc::absolute(phi, U), p)
+ fvc::ddt(rho, K) + fvc::div(rhoPhi, K)
)()() - contErr*K
)
*(
alpha1/mixture.thermo1().Cv()
+ alpha2/mixture.thermo2().Cv()
)()()
==
fvOptions(rho, T)
+ pos(Srho)
*(
surfaceFilm.Sh()()/mixture.thermo1().Cp()()
+ surfaceFilm.Tref*Srho
)
);
TEqn.relax();
fvOptions.constrain(TEqn);
TEqn.solve();
fvOptions.correct(T);
mixture.correctThermo();
mixture.correct();
}
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 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 "VoFPatchTransfer.H"
#include "twoPhaseMixtureThermo.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace regionModels
{
namespace surfaceFilmModels
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
defineTypeNameAndDebug(VoFPatchTransfer, 0);
addToRunTimeSelectionTable(transferModel, VoFPatchTransfer, dictionary);
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
VoFPatchTransfer::VoFPatchTransfer
(
surfaceFilmModel& film,
const dictionary& dict
)
:
transferModel(type(), film, dict),
deltaFactorToVoF_
(
coeffDict_.lookupOrDefault<scalar>("deltaFactorToVoF", 1.0)
),
deltaFactorToFilm_
(
coeffDict_.lookupOrDefault<scalar>("deltaFactorToFilm", 0.5)
),
alphaToVoF_
(
coeffDict_.lookupOrDefault<scalar>("alphaToVoF", 0.5)
),
alphaToFilm_
(
coeffDict_.lookupOrDefault<scalar>("alphaToFilm", 0.1)
),
transferRateCoeff_
(
coeffDict_.lookupOrDefault<scalar>("transferRateCoeff", 0.1)
)
{
const polyBoundaryMesh& pbm = film.regionMesh().boundaryMesh();
patchIDs_.setSize
(
pbm.size() - film.regionMesh().globalData().processorPatches().size()
);
if (coeffDict_.found("patches"))
{
const wordReList patchNames(coeffDict_.lookup("patches"));
const labelHashSet patchSet = pbm.patchSet(patchNames);
Info<< " applying to patches:" << nl;
label pidi = 0;
forAllConstIter(labelHashSet, patchSet, iter)
{
const label patchi = iter.key();
patchIDs_[pidi++] = patchi;
Info<< " " << pbm[patchi].name() << endl;
}
patchIDs_.setSize(pidi);
patchTransferredMasses_.setSize(pidi, 0);
}
else
{
Info<< " applying to all patches" << endl;
forAll(patchIDs_, patchi)
{
patchIDs_[patchi] = patchi;
}
patchTransferredMasses_.setSize(patchIDs_.size(), 0);
}
if (!patchIDs_.size())
{
FatalErrorInFunction
<< "No patches selected"
<< exit(FatalError);
}
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
VoFPatchTransfer::~VoFPatchTransfer()
{}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
void VoFPatchTransfer::correct
(
scalarField& availableMass,
scalarField& massToTransfer
)
{
NotImplemented;
}
void VoFPatchTransfer::correct
(
scalarField& availableMass,
scalarField& massToTransfer,
scalarField& energyToTransfer
)
{
// Do not correct if no patches selected
if (!patchIDs_.size()) return;
const scalarField& delta = film().delta();
const scalarField& rho = film().rho();
const scalarField& magSf = film().magSf();
const polyBoundaryMesh& pbm = film().regionMesh().boundaryMesh();
const twoPhaseMixtureThermo& thermo
(
film().primaryMesh().lookupObject<twoPhaseMixtureThermo>
(
twoPhaseMixtureThermo::dictName
)
);
const volScalarField& heVoF = thermo.thermo1().he();
const volScalarField& TVoF = thermo.thermo1().T();
const volScalarField CpVoF(thermo.thermo1().Cp());
const volScalarField& rhoVoF = thermo.thermo1().rho()();
const volScalarField& alphaVoF = thermo.alpha1();
forAll(patchIDs_, pidi)
{
const label patchi = patchIDs_[pidi];
label primaryPatchi = -1;
forAll(film().intCoupledPatchIDs(), i)
{
const label filmPatchi = film().intCoupledPatchIDs()[i];
if (filmPatchi == patchi)
{
primaryPatchi = film().primaryPatchIDs()[i];
}
}
if (primaryPatchi != -1)
{
scalarField deltaCoeffs
(
film().primaryMesh().boundary()[primaryPatchi].deltaCoeffs()
);
film().toRegion(patchi, deltaCoeffs);
scalarField hp(heVoF.boundaryField()[primaryPatchi]);
film().toRegion(patchi, hp);
scalarField Tp(TVoF.boundaryField()[primaryPatchi]);
film().toRegion(patchi, Tp);
scalarField Cpp(CpVoF.boundaryField()[primaryPatchi]);
film().toRegion(patchi, Cpp);
scalarField rhop(rhoVoF.boundaryField()[primaryPatchi]);
film().toRegion(patchi, rhop);
scalarField alphap(alphaVoF.boundaryField()[primaryPatchi]);
film().toRegion(patchi, alphap);
scalarField Vp
(
film().primaryMesh().boundary()[primaryPatchi]
.patchInternalField(film().primaryMesh().V())
);
film().toRegion(patchi, Vp);
const polyPatch& pp = pbm[patchi];
const labelList& faceCells = pp.faceCells();
// Accumulate the total mass removed from patch
scalar dMassPatch = 0;
forAll(faceCells, facei)
{
const label celli = faceCells[facei];
scalar dMass = 0;
if
(
delta[celli] > 2*deltaFactorToVoF_/deltaCoeffs[facei]
|| alphap[facei] > alphaToVoF_
)
{
dMass =
transferRateCoeff_*delta[celli]*rho[celli]*magSf[celli];
massToTransfer[celli] += dMass;
energyToTransfer[celli] += dMass*film().hs()[celli];
}
if
(
alphap[facei] > 0
&& delta[celli] > 0
&& delta[celli] < 2*deltaFactorToFilm_/deltaCoeffs[facei]
&& alphap[facei] < alphaToFilm_
)
{
dMass =
-transferRateCoeff_*alphap[facei]*rhop[facei]*Vp[facei];
massToTransfer[celli] += dMass;
energyToTransfer[celli] += dMass*hp[facei];
}
availableMass[celli] -= dMass;
dMassPatch += dMass;
}
patchTransferredMasses_[pidi] += dMassPatch;
addToTransferredMass(dMassPatch);
}
}
transferModel::correct();
if (writeTime())
{
scalarField patchTransferredMasses0
(
getModelProperty<scalarField>
(
"patchTransferredMasses",
scalarField(patchTransferredMasses_.size(), 0)
)
);
scalarField patchTransferredMassTotals(patchTransferredMasses_);
Pstream::listCombineGather
(
patchTransferredMassTotals,
plusEqOp<scalar>()
);
patchTransferredMasses0 += patchTransferredMassTotals;
setModelProperty<scalarField>
(
"patchTransferredMasses",
patchTransferredMasses0
);
patchTransferredMasses_ = 0;
}
}
void VoFPatchTransfer::patchTransferredMassTotals
(
scalarField& patchMasses
) const
{
// Do not correct if no patches selected
if (!patchIDs_.size()) return;
scalarField patchTransferredMasses
(
getModelProperty<scalarField>
(
"patchTransferredMasses",
scalarField(patchTransferredMasses_.size(), 0)
)
);
scalarField patchTransferredMassTotals(patchTransferredMasses_);
Pstream::listCombineGather(patchTransferredMassTotals, plusEqOp<scalar>());
forAll(patchIDs_, pidi)
{
const label patchi = patchIDs_[pidi];
patchMasses[patchi] +=
patchTransferredMasses[pidi] + patchTransferredMassTotals[pidi];
}
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace surfaceFilmModels
} // End namespace regionModels
} // End namespace Foam
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 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::regionModels::surfaceFilmModels::VoFPatchTransfer
Description
Transfer mass between the film and the VoF in the continuous phase.
SourceFiles
VoFPatchTransfer.C
\*---------------------------------------------------------------------------*/
#ifndef VoFPatchTransfer_H
#define VoFPatchTransfer_H
#include "transferModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace regionModels
{
namespace surfaceFilmModels
{
/*---------------------------------------------------------------------------*\
Class VoFPatchTransfer Declaration
\*---------------------------------------------------------------------------*/
class VoFPatchTransfer
:
public transferModel
{
// Private member functions
//- Disallow default bitwise copy construct
VoFPatchTransfer(const VoFPatchTransfer&);
//- Disallow default bitwise assignment
void operator=(const VoFPatchTransfer&);
protected:
//- Factor of the cell height above which the film is transferred
// to the VoF
scalar deltaFactorToVoF_;
//- Factor of the cell height below which the VoF may be transferred
// to the film
scalar deltaFactorToFilm_;
//- VoF limit above which all of the film is transferred to the VoF
scalar alphaToVoF_;
//- VoF limit below which the VoF may be transferred to the film
scalar alphaToFilm_;
//- Transfer rate coefficient
scalar transferRateCoeff_;
//- List of patch IDs at which the film is removed
labelList patchIDs_;
//- Transferred mass for each patch at which the film is removed
scalarField patchTransferredMasses_;
public:
//- Runtime type information
TypeName("VoFPatchTransfer");
// Constructors
//- Construct from surface film model
VoFPatchTransfer(surfaceFilmModel& film, const dictionary& dict);
//- Destructor
virtual ~VoFPatchTransfer();
// Member Functions
//- Correct
virtual void correct
(
scalarField& availableMass,
scalarField& massToTransfer
);
//- Correct kinematic and thermodynamic transfers
virtual void correct
(
scalarField& availableMass,
scalarField& massToTransfer,
scalarField& energyToTransfer
);
//- Accumulate the total mass injected for the patches into the
// scalarField provided
virtual void patchTransferredMassTotals
(
scalarField& patchMasses
) const;
};