ENH: Added new solidification and melting fvOption

This commit is contained in:
andy
2014-10-16 10:21:01 +01:00
committed by Andrew Heather
parent ade2801158
commit 5782999e7c
5 changed files with 779 additions and 1 deletions

View File

@ -12,6 +12,7 @@ $(generalSources)/semiImplicitSource/semiImplicitSource.C
derivedSources=sources/derived
$(derivedSources)/actuationDiskSource/actuationDiskSource.C
$(derivedSources)/effectivenessHeatExchangerSource/effectivenessHeatExchangerSource.C
$(derivedSources)/explicitPorositySource/explicitPorositySource.C
$(derivedSources)/MRFSource/MRFSource.C
$(derivedSources)/pressureGradientExplicitSource/pressureGradientExplicitSource.C
@ -27,7 +28,8 @@ $(derivedSources)/rotorDiskSource/trimModel/trimModel/trimModel.C
$(derivedSources)/rotorDiskSource/trimModel/trimModel/trimModelNew.C
$(derivedSources)/rotorDiskSource/trimModel/fixed/fixedTrim.C
$(derivedSources)/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.C
$(derivedSources)/effectivenessHeatExchangerSource/effectivenessHeatExchangerSource.C
$(derivedSources)/solidificationMeltingSource/solidificationMeltingSource.C
$(derivedSources)/solidificationMeltingSource/solidificationMeltingSourceIO.C
interRegion = sources/interRegion
$(interRegion)/interRegionHeatTransferModel/constantHeatTransfer/constantHeatTransfer.C

View File

@ -0,0 +1,358 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013-2014 OpenCFD Ltd.
\\/ 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 "solidificationMeltingSource.H"
#include "fvMatrices.H"
#include "basicThermo.H"
#include "uniformDimensionedFields.H"
#include "zeroGradientFvPatchFields.H"
#include "addToRunTimeSelectionTable.H"
#include "geometricOneField.H"
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
namespace Foam
{
template<>
const char* NamedEnum
<
fv::solidificationMeltingSource::thermoMode,
2
>::names[] =
{
"thermo",
"lookup"
};
namespace fv
{
defineTypeNameAndDebug(solidificationMeltingSource, 0);
addToRunTimeSelectionTable
(
option,
solidificationMeltingSource,
dictionary
);
}
}
const Foam::NamedEnum<Foam::fv::solidificationMeltingSource::thermoMode, 2>
Foam::fv::solidificationMeltingSource::thermoModeTypeNames_;
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
bool Foam::fv::solidificationMeltingSource::solveField
(
const word& fieldName
) const
{
bool result = true;
switch (mode_)
{
case mdThermo:
{
const basicThermo& thermo =
mesh_.lookupObject<basicThermo>("thermophysicalProperties");
if (fieldName != thermo.he().name())
{
result = false;
}
break;
}
case mdLookup:
{
if (fieldName != TName_)
{
result = false;
}
break;
}
default:
{
FatalErrorIn
(
"bool Foam::fv::solidificationMeltingSource::solveField"
"("
"const word&"
") const"
)
<< "Unhandled thermo mode: " << thermoModeTypeNames_[mode_]
<< abort(FatalError);
}
}
return result;
}
Foam::tmp<Foam::volScalarField>
Foam::fv::solidificationMeltingSource::Cp() const
{
switch (mode_)
{
case mdThermo:
{
const basicThermo& thermo =
mesh_.lookupObject<basicThermo>("thermophysicalProperties");
return thermo.Cp();
break;
}
case mdLookup:
{
if (CpName_ == "CpRef")
{
scalar CpRef = readScalar(coeffs_.lookup("CpRef"));
return tmp<volScalarField>
(
new volScalarField
(
IOobject
(
name_ + ":Cp",
mesh_.time().timeName(),
mesh_,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh_,
dimensionedScalar
(
"Cp",
dimEnergy/dimMass/dimTemperature,
CpRef
),
zeroGradientFvPatchScalarField::typeName
)
);
}
else
{
return mesh_.lookupObject<volScalarField>(CpName_);
}
break;
}
default:
{
FatalErrorIn
(
"Foam::tmp<Foam::volScalarField> "
"Foam::fv::solidificationMeltingSource::Cp() const"
)
<< "Unhandled thermo mode: " << thermoModeTypeNames_[mode_]
<< abort(FatalError);
}
}
return tmp<volScalarField>(NULL);
}
Foam::vector Foam::fv::solidificationMeltingSource::g() const
{
if (mesh_.foundObject<uniformDimensionedVectorField>("g"))
{
const uniformDimensionedVectorField& value =
mesh_.lookupObject<uniformDimensionedVectorField>("g");
return value.value();
}
else
{
return coeffs_.lookup("g");
}
}
void Foam::fv::solidificationMeltingSource::update(const volScalarField& Cp)
{
if (curTimeIndex_ == mesh_.time().timeIndex())
{
return;
}
if (debug)
{
Info<< type() << ": " << name_ << " - updating phase indicator" << endl;
}
// update old time alpha1 field
alpha1_.oldTime();
const volScalarField& T = mesh_.lookupObject<volScalarField>(TName_);
forAll(cells_, i)
{
label cellI = cells_[i];
scalar Tc = T[cellI];
scalar Cpc = Cp[cellI];
scalar alpha1New = alpha1_[cellI] + relax_*Cpc*(Tc - Tmelt_)/L_;
alpha1_[cellI] = max(0, min(alpha1New, 1));
deltaT_[i] = Tc - Tmelt_;
}
alpha1_.correctBoundaryConditions();
curTimeIndex_ = mesh_.time().timeIndex();
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::fv::solidificationMeltingSource::solidificationMeltingSource
(
const word& sourceName,
const word& modelType,
const dictionary& dict,
const fvMesh& mesh
)
:
option(sourceName, modelType, dict, mesh),
Tmelt_(readScalar(coeffs_.lookup("Tmelt"))),
L_(readScalar(coeffs_.lookup("L"))),
relax_(coeffs_.lookupOrDefault("relax", 0.9)),
mode_(thermoModeTypeNames_.read(coeffs_.lookup("thermoMode"))),
rhoRef_(readScalar(coeffs_.lookup("rhoRef"))),
TName_(coeffs_.lookupOrDefault<word>("TName", "T")),
CpName_(coeffs_.lookupOrDefault<word>("CpName", "Cp")),
UName_(coeffs_.lookupOrDefault<word>("UName", "U")),
phiName_(coeffs_.lookupOrDefault<word>("phiName", "phi")),
Cu_(coeffs_.lookupOrDefault<scalar>("Cu", 100000)),
q_(coeffs_.lookupOrDefault("q", 0.001)),
beta_(readScalar(coeffs_.lookup("beta"))),
alpha1_
(
IOobject
(
name_ + ":alpha1",
mesh.time().timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("alpha1", dimless, 0.0),
zeroGradientFvPatchScalarField::typeName
),
curTimeIndex_(-1),
deltaT_(cells_.size(), 0)
{
fieldNames_.setSize(1, "source");
applied_.setSize(1, false);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::fv::solidificationMeltingSource::alwaysApply() const
{
return true;
}
void Foam::fv::solidificationMeltingSource::addSup
(
fvMatrix<scalar>& eqn,
const label fieldI
)
{
apply(geometricOneField(), eqn);
}
void Foam::fv::solidificationMeltingSource::addSup
(
const volScalarField& rho,
fvMatrix<scalar>& eqn,
const label fieldI
)
{
apply(rho, eqn);
}
void Foam::fv::solidificationMeltingSource::addSup
(
fvMatrix<vector>& eqn,
const label fieldI
)
{
const volVectorField& U = eqn.psi();
if (U.name() != UName_)
{
return;
}
if (debug)
{
Info<< type() << ": applying source to " << UName_ << endl;
}
const volScalarField Cp(this->Cp());
update(Cp);
vector g = this->g();
scalarField& Sp = eqn.diag();
vectorField& Su = eqn.source();
const scalarField& V = mesh_.V();
forAll(cells_, i)
{
label cellI = cells_[i];
scalar Vc = V[cellI];
scalar alpha1c = alpha1_[cellI];
scalar S = -Cu_*sqr(1.0 - alpha1c)/(pow3(alpha1c) + q_);
vector Sb = rhoRef_*g*beta_*deltaT_[i];
Sp[cellI] += Vc*S;
Su[cellI] += Vc*Sb;
}
}
void Foam::fv::solidificationMeltingSource::addSup
(
const volScalarField& rho,
fvMatrix<vector>& eqn,
const label fieldI
)
{
// momentum source uses a Boussinesq approximation - redirect
addSup(eqn, fieldI);
}
// ************************************************************************* //

View File

@ -0,0 +1,280 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013-2014 OpenCFD Ltd.
\\/ 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::fv::solidificationMeltingSource
Description
This source is designed to model the effect of solidification and melting
processes, e.g. windhield defrosting. The phase change occurs at the
melting temperature, \c Tmelt.
The presence of the solid phase in the flow field is incorporated into the
model as a momentum porosity contribution; the energy associated with the
phase change is added as an enthalpy contribution.
Based on the references:
1. V.R. Voller and C. Prakash, A fixed grid numerical modelling
methodology for convection-diffusion mushy phase-change problems,
Int. J. Heat Mass Transfer 30(8):17091719, 1987.
2. C.R. Swaminathan. and V.R. Voller, A general enthalpy model for
modeling solidification processes, Metallurgical Transactions
23B:651664, 1992.
The model generates a field \c \<name\>:alpha1 which can be visualised to
to show the melt distribution as a fraction [0-1]
\heading Source usage
Example usage:
\verbatim
solidificationMeltingSource1
{
type solidificationMeltingSource;
active on;
selectionMode cellZone;
cellZone iceZone;
solidificationMeltingSourceCoeffs
{
Tmelt 273;
L 334000;
thermoMode thermo;
beta 50e-6;
rhoRef 800;
}
}
\endverbatim
Where:
\table
Property | Description | Required | Default value
Tmelt | Melting temperature [K] | yes |
L | Latent heat of fusion [J/kg] | yes |
relax | Relaxation coefficient [0-1] | no | 0.9
thermoMode | Thermo mode [thermo|lookup] | yes |
rhoRef | Reference (solid) density | yes |
rhoName | Name of density field | no | rho
TName | Name of temperature field | no | T
CpName | Name of specific heat capacity field | no | Cp
UName | Name of velocity field | no | U
phiName | Name of flux field | no | phi
Cu | Model coefficient | no | 100000
q | Model coefficient | no | 0.001
beta | Thermal expansion coefficient [1/K] | yes |
g | Accelerartion due to gravity | no |
\endtable
SourceFiles
solidificationMeltingSource.C
solidificationMeltingSourceIO.C
\*---------------------------------------------------------------------------*/
#ifndef solidificationMeltingSource_H
#define solidificationMeltingSource_H
#include "fvMesh.H"
#include "volFields.H"
#include "fvOption.H"
#include "NamedEnum.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace fv
{
/*---------------------------------------------------------------------------*\
Class solidificationMeltingSource Declaration
\*---------------------------------------------------------------------------*/
class solidificationMeltingSource
:
public option
{
public:
enum thermoMode
{
mdThermo,
mdLookup
};
static const NamedEnum<thermoMode, 2> thermoModeTypeNames_;
private:
// Private data
//- Temperature at which melting occurs [K]
scalar Tmelt_;
//- Latent heat of fusion [J/kg]
scalar L_;
//- Phase fraction under-relaxation coefficient
scalar relax_;
//- Thermodynamics mode
thermoMode mode_;
//- Reference density - typically the solid density
scalar rhoRef_;
//- Name of temperature field - default = "T" (optional)
word TName_;
//- Name of specific heat capacity field - default = "Cp" (optional)
word CpName_;
//- Name of velocity field - default = "U" (optional)
word UName_;
//- Name of flux field - default = "phi" (optional)
word phiName_;
//- Mushy region momentum sink coefficient [1/s]; default = 10^5
scalar Cu_;
//- Coefficient used in porosity calc - default = 0.001
scalar q_;
//- Thermal expansion coefficient [1/K]
scalar beta_;
//- Phase fraction indicator field
volScalarField alpha1_;
//- Current time index (used for updating)
label curTimeIndex_;
//- Temperature change cached for source calculation when alpha1 updated
scalarField deltaT_;
// Private Member Functions
//- Flag to indicate whether to solve for given field
bool solveField(const word& fieldName) const;
//- Return the specific heat capacity field
tmp<volScalarField> Cp() const;
//- Return the gravity vector
vector g() const;
//- Update the model
void update(const volScalarField& Cp);
//- Helper function to apply to the energy equation
template<class RhoFieldType>
void apply(const RhoFieldType& rho, fvMatrix<scalar>& eqn);
//- Disallow default bitwise copy construct
solidificationMeltingSource(const solidificationMeltingSource&);
//- Disallow default bitwise assignment
void operator=(const solidificationMeltingSource&);
public:
//- Runtime type information
TypeName("solidificationMeltingSource");
// Constructors
//- Construct from explicit source name and mesh
solidificationMeltingSource
(
const word& sourceName,
const word& modelType,
const dictionary& dict,
const fvMesh& mesh
);
// Member Functions
//- Flag to bypass the apply flag list checking
virtual bool alwaysApply() const;
// Add explicit and implicit contributions
//- Add explicit contribution to enthalpy equation
virtual void addSup(fvMatrix<scalar>& eqn, const label fieldI);
//- Add implicit contribution to momentum equation
virtual void addSup(fvMatrix<vector>& eqn, const label fieldI);
// Add explicit and implicit contributions to compressible equation
//- Add explicit contribution to compressible enthalpy equation
virtual void addSup
(
const volScalarField& rho,
fvMatrix<scalar>& eqn,
const label fieldI
);
//- Add implicit contribution to compressible momentum equation
virtual void addSup
(
const volScalarField& rho,
fvMatrix<vector>& eqn,
const label fieldI
);
// I-O
//- Write the source properties
virtual void writeData(Ostream&) const;
//- Read source dictionary
virtual bool read(const dictionary& dict);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace fv
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "solidificationMeltingSourceTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,68 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013-2014 OpenCFD Ltd.
\\/ 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 "solidificationMeltingSource.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::fv::solidificationMeltingSource::writeData(Ostream& os) const
{
os << indent << name_ << endl;
dict_.write(os);
}
bool Foam::fv::solidificationMeltingSource::read(const dictionary& dict)
{
if (option::read(dict))
{
coeffs_.lookup("Tmelt") >> Tmelt_;
coeffs_.lookup("L") >> L_;
coeffs_.readIfPresent("relax", relax_);
mode_ = thermoModeTypeNames_.read(coeffs_.lookup("thermoMode"));
coeffs_.lookup("rhoRef") >> rhoRef_;
coeffs_.readIfPresent("TName", TName_);
coeffs_.readIfPresent("UName", UName_);
coeffs_.readIfPresent("Cu", Cu_);
coeffs_.readIfPresent("q", q_);
coeffs_.readIfPresent("beta", beta_);
return true;
}
else
{
return false;
}
return false;
}
// ************************************************************************* //

View File

@ -0,0 +1,70 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2014 OpenCFD Ltd.
\\/ 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 "fvMatrices.H"
#include "fvcDdt.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class RhoFieldType>
void Foam::fv::solidificationMeltingSource::apply
(
const RhoFieldType& rho,
fvMatrix<scalar>& eqn
)
{
if (!solveField(eqn.psi().name()))
{
return;
}
if (debug)
{
Info<< type() << ": applying source to " << eqn.psi().name() << endl;
}
const volScalarField Cp(this->Cp());
update(Cp);
dimensionedScalar L("L", dimEnergy/dimMass, L_);
// contributions added to rhs of solver equation
if (eqn.psi().dimensions() == dimTemperature)
{
// isothermal phase change - only include time derivative
// eqn -= L/Cp*(fvc::ddt(rho, alpha1_) + fvc::div(phi, alpha1_));
eqn -= L/Cp*(fvc::ddt(rho, alpha1_));
}
else
{
// isothermal phase change - only include time derivative
// eqn -= L*(fvc::ddt(rho, alpha1_) + fvc::div(phi, alpha1_));
eqn -= L*(fvc::ddt(rho, alpha1_));
}
}
// ************************************************************************* //