phaseHydrostaticPressureFvPatchScalarField.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 | Copyright (C) 2011-2017 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
28 #include "fvPatchFieldMapper.H"
29 #include "volFields.H"
30 #include "surfaceFields.H"
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
37 (
38  const fvPatch& p,
40 )
41 :
42  mixedFvPatchScalarField(p, iF),
43  phaseFraction_("alpha"),
44  rho_(0.0),
45  pRefValue_(0.0),
46  pRefPoint_(Zero)
47 {
48  this->refValue() = 0.0;
49  this->refGrad() = 0.0;
50  this->valueFraction() = 0.0;
51 }
52 
53 
56 (
57  const fvPatch& p,
59  const dictionary& dict
60 )
61 :
62  mixedFvPatchScalarField(p, iF),
63  phaseFraction_(dict.lookupOrDefault<word>("phaseFraction", "alpha")),
64  rho_(readScalar(dict.lookup("rho"))),
65  pRefValue_(readScalar(dict.lookup("pRefValue"))),
66  pRefPoint_(dict.lookup("pRefPoint"))
67 {
68  this->refValue() = pRefValue_;
69 
70  if (dict.found("value"))
71  {
72  fvPatchScalarField::operator=
73  (
74  scalarField("value", dict, p.size())
75  );
76  }
77  else
78  {
79  fvPatchScalarField::operator=(this->refValue());
80  }
81 
82  this->refGrad() = 0.0;
83  this->valueFraction() = 0.0;
84 }
85 
86 
89 (
91  const fvPatch& p,
93  const fvPatchFieldMapper& mapper
94 )
95 :
96  mixedFvPatchScalarField(ptf, p, iF, mapper),
97  phaseFraction_(ptf.phaseFraction_),
98  rho_(ptf.rho_),
99  pRefValue_(ptf.pRefValue_),
100  pRefPoint_(ptf.pRefPoint_)
101 {}
102 
103 
106 (
108 )
109 :
110  mixedFvPatchScalarField(ptf),
111  phaseFraction_(ptf.phaseFraction_)
112 {}
113 
114 
117 (
120 )
121 :
122  mixedFvPatchScalarField(ptf, iF),
123  phaseFraction_(ptf.phaseFraction_),
124  rho_(ptf.rho_),
125  pRefValue_(ptf.pRefValue_),
126  pRefPoint_(ptf.pRefPoint_)
127 {}
128 
129 
130 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
131 
133 {
134  if (this->updated())
135  {
136  return;
137  }
138 
139  const scalarField& alphap =
140  patch().lookupPatchField<volScalarField, scalar>
141  (
143  );
144 
146  db().lookupObject<uniformDimensionedVectorField>("g");
147 
148  // scalar rhor = 1000;
149  // scalarField alphap1 = max(min(alphap, 1.0), 0.0);
150  // valueFraction() = alphap1/(alphap1 + rhor*(1.0 - alphap1));
151  valueFraction() = max(min(alphap, scalar(1)), scalar(0));
152 
153  refValue() =
154  pRefValue_
155  + rho_*((g.value() & patch().Cf()) - (g.value() & pRefPoint_));
156 
157  mixedFvPatchScalarField::updateCoeffs();
158 }
159 
160 
162 {
164  if (phaseFraction_ != "alpha")
165  {
166  os.writeKeyword("phaseFraction")
168  }
169  os.writeKeyword("rho") << rho_ << token::END_STATEMENT << nl;
170  os.writeKeyword("pRefValue") << pRefValue_ << token::END_STATEMENT << nl;
171  os.writeKeyword("pRefPoint") << pRefPoint_ << token::END_STATEMENT << nl;
172  writeEntry("value", os);
173 }
174 
175 
176 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
177 
178 void Foam::phaseHydrostaticPressureFvPatchScalarField::operator=
179 (
180  const fvPatchScalarField& ptf
181 )
182 {
184  (
185  valueFraction()*refValue() + (1 - valueFraction())*ptf
186  );
187 }
188 
189 
190 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
191 
192 namespace Foam
193 {
195  (
198  );
199 }
200 
201 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
Foam::surfaceFields.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:431
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
Macros for easy insertion into run-time selection tables.
This boundary condition provides a phase-based hydrostatic pressure condition, calculated as: ...
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
phaseHydrostaticPressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
static const zero Zero
Definition: zero.H:91
virtual label size() const
Return size.
Definition: fvPatch.H:161
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
volScalarField scalarField(fieldObject, mesh)
static const char nl
Definition: Ostream.H:262
const dimensionedVector & g
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576