filmFlux.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) 2019 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "filmFlux.H"
29#include "volFields.H"
30#include "surfaceFields.H"
31#include "surfaceInterpolate.H"
33
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
38namespace functionObjects
39{
42}
43}
44
45
47Foam::functionObjects::filmFlux::filmModel()
48{
49 return time_.objectRegistry::lookupObject<filmType>(filmName_);
50}
51
52
53// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
54
56(
57 const word& name,
58 const Time& runTime,
59 const dictionary& dict
60)
61:
63 filmName_("surfaceFilmProperties"),
64 resultName_(scopedName("filmFlux"))
65{
66 read(dict);
67}
68
69
70// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
71
73{
75 {
76 dict.readIfPresent<word>("film", filmName_);
77 dict.readIfPresent<word>("result", resultName_);
78
79 return true;
80 }
81
82 return false;
83}
84
85
87{
88 const auto& model = filmModel();
89
90 const fvMesh& filmMesh = model.regionMesh();
91
92 auto* resultPtr = filmMesh.getObjectPtr<surfaceScalarField>(resultName_);
93
94 if (!resultPtr)
95 {
96 resultPtr = new surfaceScalarField
97 (
99 (
100 resultName_,
101 time_.timeName(),
102 filmMesh,
105 ),
106 filmMesh,
108 );
109
110 filmMesh.objectRegistry::store(resultPtr);
111 }
112
113 auto& result = *resultPtr;
114
115 const surfaceScalarField& phi = model.phi();
116 const volScalarField& magSf = model.magSf();
117 const volScalarField::Internal& vol = filmMesh.V();
118
119 volScalarField height
120 (
122 (
123 "height",
124 time_.timeName(),
125 filmMesh,
128 false
129 ),
130 filmMesh,
132 zeroGradientFvPatchScalarField::typeName
133 );
134
135 auto& heightc = height.ref();
136
137 heightc = max(dimensionedScalar("eps", dimLength, ROOTVSMALL), vol/magSf());
139
140 result = phi/fvc::interpolate(height);
141
142 return true;
143}
144
145
147{
148 const auto& filmMesh = filmModel().regionMesh();
149
150 filmMesh.lookupObject<surfaceScalarField>(resultName_).write();
151
152 return true;
153}
154
155
156// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
surfaceScalarField & phi
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
virtual bool read()
Re-read model coefficients if they have changed.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Abstract base-class for Time/database function objects.
Calculates the mass flux for surface film models.
Definition: filmFlux.H:146
virtual bool execute()
Execute.
Definition: filmFlux.C:86
virtual bool write()
Write the field.
Definition: filmFlux.C:146
virtual bool read(const dictionary &)
Read the field data.
Definition: filmFlux.C:72
Base class for function objects, adding functionality to read/write state information (data required ...
const Time & time_
Reference to the time database.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
Type * getObjectPtr(const word &name, const bool recursive=false) const
Kinematic form of single-cell layer surface film model.
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
engineTime & runTime
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
runTime write()
dictionary dict
Foam::surfaceFields.