Commit 8016af37 authored by Henry Weller's avatar Henry Weller
Browse files

reactingTwoPhaseEulerFoam: Enhanced support for wall boiling

Contributed by Juho Peltola, VTT

Notable changes:

    1. The same wall function is now used for both phases, but user must
       specify phaseType ‘liquid’ or ‘vapor’

    2. Runtime selectable submodels for:
       - wall heat flux partitioning between the phases
       - nucleation site density
       - bubble departure frequency
       - bubble departure diameter

    3. An additional iteration loop for the wall boiling model in case
       the initial guess for the wall temperature proves to be poor.

The wallBoiling tutorial has been updated to demonstrate this new functionality.
parent 6c3c2f13
......@@ -36,10 +36,31 @@ 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/alphatPhaseChangeJayatillekeWallFunction/alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField.C
derivedFvPatchFields/alphatFixedDmdtWallBoilingWallFunction/alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField.C
derivedFvPatchFields/alphatWallBoilingWallFunction/alphatWallBoilingWallFunctionFvPatchScalarField.C
derivedFvPatchFields/copiedFixedValue/copiedFixedValueFvPatchScalarField.C
derivedFvPatchFields/fixedMultiPhaseHeatFlux/fixedMultiPhaseHeatFluxFvPatchScalarField.C
LIB = $(FOAM_LIBBIN)/libtwoPhaseReactingTurbulenceModels
......@@ -28,11 +28,17 @@ Group
grpCmpWallFunctions
Description
A thermal wall function for simulation of subcooled nucleate wall boiling.
A thermal wall function for simulation of subcooled nucleate wall boiling
with runtime selctable submodels for:
- wall heat flux partitioning model
- nucleation site density
- bubble departure frequency
- bubble departure diameter
Implements a version of the well-known RPI wall boiling model
(Kurul & Podowski, 1991). The model implementation is similar to the model
described by Peltola & Pättikangas (2012).
described by Peltola & Pättikangas (2012) but has been extended with the
wall heat flux partitioning models.
References:
\verbatim
......@@ -51,8 +57,66 @@ Description
Daejeon, Korea, September 10-12 2012
\endverbatim
See also
Foam::fixedValueFvPatchField
Usage
\table
Property | Description | Required | Default value
phaseType | 'vapor' or 'liquid' | yes |
relax |wall boiling model relaxation| yes |
Prt | inherited from alphatPhaseChangeJayatillekeWallFunction
Cmu | inherited from alphatPhaseChangeJayatillekeWallFunction
kappa | inherited from alphatPhaseChangeJayatillekeWallFunction
E | inherited from alphatPhaseChangeJayatillekeWallFunction
dmdt | phase change mass flux | yes |
value | initial alphat value | yes |
if phaseType 'vapor':
partitioningModel| | yes |
if phaseType 'liquid':
partitioningModel| | yes |
nucleationSiteModel| | yes |
departureDiamModel| | yes |
departureFreqModel| | yes |
\endtable
NOTE: Runtime selectabale submodels may require model specific entries
Example usage:
\verbatim
hotWall
{
type compressible::alphatWallBoiling2WallFunction;
phaseType liquid;
Prt 0.85;
Cmu 0.09;
kappa 0.41;
E 9.8;
relax 0.001;
dmdt uniform 0;
partitioningModel
{
type Lavieville;
alphaCrit 0.2;
}
nucleationSiteModel
{
type LemmertChawla;
}
departureDiamModel
{
type TolubinskiKostanchuk;
}
departureFreqModel
{
type Cole;
}
value uniform 0.01;
\endverbatim
SeeAlso
Foam::alphatPhaseChangeJayatillekeWallFunctionFvPatchField
SourceFiles
alphatWallBoilingWallFunctionFvPatchScalarField.C
......@@ -63,6 +127,10 @@ SourceFiles
#define compressiblealphatWallBoilingWallFunctionFvPatchScalarField_H
#include "alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField.H"
#include "partitioningModel.H"
#include "nucleationSiteModel.H"
#include "departureDiameterModel.H"
#include "departureFrequencyModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......@@ -79,8 +147,28 @@ class alphatWallBoilingWallFunctionFvPatchScalarField
:
public alphatPhaseChangeJayatillekeWallFunctionFvPatchScalarField
{
public:
// Data types
//- Enumeration listing the possible operational modes
enum phaseType
{
vaporPhase,
liquidPhase
};
private:
// Private data
//- Heat source type names
static const NamedEnum<phaseType, 2> phaseTypeNames_;
//- Heat source type
phaseType phaseType_;
//- dmdt relaxationFactor
scalar relax_;
......@@ -90,6 +178,22 @@ class alphatWallBoilingWallFunctionFvPatchScalarField
//- Convective turbulent thermal diffusivity
scalarField alphatConv_;
//- Run-time selected heat flux partitioning model
autoPtr<wallBoilingModels::partitioningModel>
partitioningModel_;
//- Run-time selected nucleation site density model
autoPtr<wallBoilingModels::nucleationSiteModel>
nucleationSiteModel_;
//- Run-time selected bubble departure diameter model
autoPtr<wallBoilingModels::departureDiameterModel>
departureDiamModel_;
//- Run-time selected bubble departure frequency model
autoPtr<wallBoilingModels::departureFrequencyModel>
departureFreqModel_;
public:
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 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 "KocamustafaogullariIshii.H"
#include "addToRunTimeSelectionTable.H"
#include "uniformDimensionedFields.H"
#include "compressibleTurbulenceModel.H"
#include "ThermalDiffusivity.H"
#include "PhaseCompressibleTurbulenceModel.H"
#include "phaseSystem.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace wallBoilingModels
{
namespace departureDiameterModels
{
defineTypeNameAndDebug(KocamustafaogullariIshii, 0);
addToRunTimeSelectionTable
(
departureDiameterModel,
KocamustafaogullariIshii,
dictionary
);
}
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::wallBoilingModels::departureDiameterModels::
KocamustafaogullariIshii::KocamustafaogullariIshii
(
const dictionary& dict
)
:
departureDiameterModel(),
phi_(readScalar(dict.lookup("phi")))
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::wallBoilingModels::departureDiameterModels::
KocamustafaogullariIshii::~KocamustafaogullariIshii()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::tmp<Foam::scalarField>
Foam::wallBoilingModels::departureDiameterModels::
KocamustafaogullariIshii::dDeparture
(
const phaseModel& liquid,
const phaseModel& vapor,
const label patchi,
const scalarField& Tsub
) const
{
// Gravitational acceleration
const uniformDimensionedVectorField& g =
liquid.mesh().lookupObject<uniformDimensionedVectorField>("g");
const fvPatchScalarField& rhoLiquid =
liquid.turbulence().rho().boundaryField()[patchi];
const fvPatchScalarField& rhoVapor =
vapor.turbulence().rho().boundaryField()[patchi];
const scalarField rhoM((rhoLiquid-rhoVapor)/rhoVapor);
const tmp<volScalarField>& tsigma
(
liquid.fluid().sigma(phasePairKey(liquid.name(),vapor.name()))
);
const volScalarField& sigma = tsigma();
const fvPatchScalarField& sigmaw = sigma.boundaryField()[patchi];
return
0.0012*pow(rhoM,0.9)*0.0208*phi_
*sqrt(sigmaw/(mag(g.value())*(rhoLiquid-rhoVapor)));
}
void Foam::wallBoilingModels::departureDiameterModels::
KocamustafaogullariIshii::write(Ostream& os) const
{
departureDiameterModel::write(os);
os.writeKeyword("phi") << phi_ << token::END_STATEMENT << nl;
}
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 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::wallBoilingModels::departureDiameterModels::KocamustafaogullariIshii
Description
A correlation for bubble departure diameter.
Requires model parameter 'phi': contact angle in degrees.
Reference:
\verbatim
Kocamustafaogullari, G., & Ishii, M. (1983).
Interfacial area and nucleation site density in boiling systems.
International Journal of Heat and Mass Transfer, 26(9), 1377-1387.
\endverbatim
SourceFiles
KocamustafaogullariIshii.C
\*---------------------------------------------------------------------------*/
#ifndef KocamustafaogullariIshii_H
#define KocamustafaogullariIshii_H
#include "departureDiameterModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace wallBoilingModels
{
namespace departureDiameterModels
{
/*---------------------------------------------------------------------------*\
Class KocamustafaogullariIshii Declaration
\*---------------------------------------------------------------------------*/
class KocamustafaogullariIshii
:
public departureDiameterModel
{
// Private data
//- Contact angle
scalar phi_;
public:
//- Runtime type information
TypeName("KocamustafaogullariIshii");
// Constructors
//- Construct from a dictionary
KocamustafaogullariIshii(const dictionary& dict);
//- Destructor
virtual ~KocamustafaogullariIshii();
// Member Functions
//- Calculate and return the departure diameter field
virtual tmp<scalarField> dDeparture
(
const phaseModel& liquid,
const phaseModel& vapor,
const label patchi,
const scalarField& Tsub
) const;
virtual void write(Ostream& os) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace departureDiameterModels
} // End namespace wallBoilingModels
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 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 "TolubinskiKostanchuk.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace wallBoilingModels
{
namespace departureDiameterModels
{
defineTypeNameAndDebug(TolubinskiKostanchuk, 0);
addToRunTimeSelectionTable
(
departureDiameterModel,
TolubinskiKostanchuk,
dictionary
);
}
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::wallBoilingModels::departureDiameterModels::
TolubinskiKostanchuk::TolubinskiKostanchuk
(
const dictionary& dict
)
:
departureDiameterModel()
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::wallBoilingModels::departureDiameterModels::
TolubinskiKostanchuk::~TolubinskiKostanchuk()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::tmp<Foam::scalarField>
Foam::wallBoilingModels::departureDiameterModels::
TolubinskiKostanchuk::dDeparture
(
const phaseModel& liquid,
const phaseModel& vapor,
const label patchi,
const scalarField& Tsub
) const
{
return max(min(0.0006*exp(-Tsub/45), scalar(0.0014)), scalar(1e-6));
}
// ************************************************************************* //