liquidPropertiesSurfaceTension.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017 OpenFOAM Foundation
9 Copyright (C) 2019 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
30#include "liquidThermo.H"
32
33// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34
35namespace Foam
36{
37namespace surfaceTensionModels
38{
39 defineTypeNameAndDebug(liquidProperties, 0);
41 (
42 surfaceTensionModel,
43 liquidProperties,
44 dictionary
45 );
46}
47}
48
49
50// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51
53(
54 const dictionary& dict,
55 const fvMesh& mesh
56)
57:
58 surfaceTensionModel(mesh),
59 phaseName_(dict.get<word>("phase"))
60{}
61
62
63// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
64
67{
70 (
72 );
73
74 const Foam::liquidProperties& liquid = thermo.mixture().properties();
75
76 tmp<volScalarField> tsigma
77 (
79 (
80 IOobject
81 (
82 "sigma",
83 mesh_.time().timeName(),
84 mesh_,
87 false
88 ),
89 mesh_,
90 dimSigma
91 )
92 );
93 volScalarField& sigma = tsigma.ref();
94
95 const volScalarField& T = thermo.T();
96 const volScalarField& p = thermo.p();
97
100 const volScalarField::Internal& Ti = T;
101
102 forAll(sigmai, celli)
103 {
104 sigmai[celli] = liquid.sigma(pi[celli], Ti[celli]);
105 }
106
107 volScalarField::Boundary& sigmaBf = sigma.boundaryFieldRef();
110
111 forAll(sigmaBf, patchi)
112 {
113 scalarField& sigmaPf = sigmaBf[patchi];
114 const scalarField& pPf = pBf[patchi];
115 const scalarField& TPf = TBf[patchi];
116
117 forAll(sigmaPf, facei)
118 {
119 sigmaPf[facei] = liquid.sigma(pPf[facei], TPf[facei]);
120 }
121 }
122
123 return tsigma;
124}
125
126
128(
129 const dictionary& dict
130)
131{
132 return true;
133}
134
135
137(
138 Ostream& os
139) const
140{
141 if (surfaceTensionModel::writeData(os))
142 {
143 return os.good();
144 }
145
146 return false;
147}
148
149
150// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
DimensionedField< scalar, volMesh > Internal
The internal field type from which this GeometricField is derived.
GeometricBoundaryField< scalar, fvPatchField, volMesh > Boundary
Type of boundary fields.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
static const word dictName
Definition: basicThermo.H:256
The thermophysical properties of a liquid.
virtual bool writeData(Ostream &os) const
Write in dictionary format.
virtual bool readDict(const dictionary &dict)
Update surface tension coefficient from given dictionary.
virtual tmp< volScalarField > sigma() const
Surface tension coefficient.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
volScalarField & p
const volScalarField & T
dynamicFvMesh & mesh
OBJstream os(runTime.globalPath()/outputName)
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
constexpr scalar pi(M_PI)
Namespace for OpenFOAM.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
heRhoThermo< rhoThermo, pureMixture< species::thermo< thermophysicalPropertiesSelector< liquidProperties >, sensibleInternalEnergy > > > heRhoThermopureMixtureliquidProperties
Definition: liquidThermo.H:55
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333