waveAlphaFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2017-2018 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 
30 #include "levelSet.H"
31 #include "surfaceFields.H"
32 #include "volFields.H"
33 #include "fvMeshSubset.H"
34 
35 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36 
38 (
39  const fvPatch& p,
41 )
42 :
43  mixedFvPatchScalarField(p, iF),
44  UName_("U"),
45  liquid_(true),
46  inletOutlet_(true)
47 {
48  refValue() = Zero;
49  refGrad() = Zero;
50  valueFraction() = 0;
51 }
52 
53 
55 (
56  const fvPatch& p,
58  const dictionary& dict
59 )
60 :
61  mixedFvPatchScalarField(p, iF),
62  UName_(dict.lookupOrDefault<word>("U", "U")),
63  liquid_(dict.lookupOrDefault<Switch>("liquid", true)),
64  inletOutlet_(dict.lookupOrDefault<Switch>("inletOutlet", true))
65 {
66  if (dict.found("value"))
67  {
68  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
69  }
70  else
71  {
72  fvPatchScalarField::operator=(patchInternalField());
73  }
74 
75  refValue() = *this;
76  refGrad() = Zero;
77  valueFraction() = 0;
78 }
79 
80 
82 (
83  const waveAlphaFvPatchScalarField& ptf,
84  const fvPatch& p,
86  const fvPatchFieldMapper& mapper
87 )
88 :
89  mixedFvPatchScalarField(ptf, p, iF, mapper),
90  UName_(ptf.UName_),
91  liquid_(ptf.liquid_),
92  inletOutlet_(ptf.inletOutlet_)
93 {}
94 
95 
97 (
99 )
100 :
101  mixedFvPatchScalarField(ptf),
102  UName_(ptf.UName_),
103  liquid_(ptf.liquid_),
104  inletOutlet_(ptf.inletOutlet_)
105 {}
106 
108 (
109  const waveAlphaFvPatchScalarField& ptf,
111 )
112 :
113  mixedFvPatchScalarField(ptf, iF),
114  UName_(ptf.UName_),
115  liquid_(ptf.liquid_),
116  inletOutlet_(ptf.inletOutlet_)
117 {}
118 
119 
120 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
121 
123 {
124  const scalar t = db().time().timeOutputValue();
125  const waveSuperposition& waves = waveSuperposition::New(db());
126 
127  return
129  (
130  patch(),
131  waves.height(t, patch().Cf()),
132  waves.height(t, patch().patch().localPoints()),
133  !liquid_
134  );
135 }
136 
137 
139 {
140  const scalar t = db().time().timeOutputValue();
141  const waveSuperposition& waves = waveSuperposition::New(db());
143  refCast<const waveVelocityFvPatchVectorField>
144  (
145  patch().lookupPatchField<volVectorField, scalar>(UName_)
146  );
147 
148  const fvMeshSubset& subset = Up.faceCellSubset();
149  const fvMesh& meshs = subset.subMesh();
150  const label patchis = findIndex(subset.patchMap(), patch().index());
151 
152  const scalarField alphas
153  (
155  (
156  meshs,
157  waves.height(t, meshs.cellCentres())(),
158  waves.height(t, meshs.points())(),
159  !liquid_
160  )
161  );
162 
163  tmp<scalarField> tResult(new scalarField(patch().size()));
164  scalarField& result = tResult.ref();
165 
166  if (patchis != -1)
167  {
168  forAll(meshs.boundary()[patchis], is)
169  {
170  const label fs = is + meshs.boundary()[patchis].patch().start();
171  const label cs = meshs.boundary()[patchis].faceCells()[is];
172  const label f = subset.faceMap()[fs];
173  const label i = patch().patch().whichFace(f);
174  result[i] = alphas[cs];
175  }
176  }
177 
178  return tResult;
179 }
180 
181 
183 {
184  if (updated())
185  {
186  return;
187  }
188 
189  const fvPatchVectorField& Up =
190  patch().lookupPatchField<volVectorField, scalar>(UName_);
191 
192  if (!isA<waveVelocityFvPatchVectorField>(Up))
193  {
195  << "The corresponding condition for the velocity "
196  << "field " << UName_ << " on patch " << patch().name()
197  << " is not of type " << waveVelocityFvPatchVectorField::typeName
198  << exit(FatalError);
199  }
200 
201  const waveVelocityFvPatchVectorField& Uwp =
202  refCast<const waveVelocityFvPatchVectorField>(Up);
203 
204  const fvPatchScalarField& pp =
205  patch().lookupPatchField<volScalarField, scalar>(Uwp.pName());
206 
207  if (isA<wavePressureFvPatchScalarField>(pp))
208  {
209  const scalarField alpha(this->alpha()), alphan(this->alphan());
210  const scalarField out(pos0(Uwp.U() & patch().Sf()));
211 
212  valueFraction() = out;
213  refValue() = alpha;
214  refGrad() = (alpha - alphan)*patch().deltaCoeffs();
215  }
216  else
217  {
218  refValue() = alpha();
219 
220  if (inletOutlet_)
221  {
222  const scalarField& phip =
223  patch().lookupPatchField<surfaceScalarField, scalar>("phi");
224  const scalarField out(pos0(phip));
225 
226  valueFraction() = 1 - out;
227  }
228  else
229  {
230  valueFraction() = 1;
231  }
232  }
233 
234  mixedFvPatchScalarField::updateCoeffs();
235 }
236 
237 
239 (
240  Ostream& os
241 ) const
242 {
244  writeEntryIfDifferent<word>(os, "U", "U", UName_);
245  writeEntryIfDifferent<Switch>(os, "inletOutlet", true, inletOutlet_);
246  writeEntryIfDifferent<Switch>(os, "liquid", true, liquid_);
247 }
248 
249 
250 // * * * * * * * * * * * * * * Build Macro Function * * * * * * * * * * * * //
251 
252 namespace Foam
253 {
255 }
256 
257 // ************************************************************************* //
Foam::surfaceFields.
A wrapper around a list of wave models. Superimposes the modelled values of elevation and velocity...
const labelList & patchMap() const
Return patch map.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const fvMeshSubset & faceCellSubset() const
Access the face-cell subset.
tmp< scalarField > alpha() const
Return the current modelled phase fraction field.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
tmp< scalarField > levelSetFraction(const fvMesh &mesh, const scalarField &levelC, const scalarField &levelP, const bool above)
Calculate the volume-fraction that a level set occupies. This gives the.
Definition: levelSet.C:34
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
This boundary condition provides a waveAlpha condition. This sets the phase fraction to that specifie...
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
Macros for easy insertion into run-time selection tables.
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
virtual tmp< scalarField > height(const scalar t, const vectorField &p) const
Get the height above the waves at a given time and global positions.
const word & pName() const
Access the name of the pressure field.
A class for handling words, derived from string.
Definition: word.H:59
virtual void write(Ostream &) const
Write.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Foam::fvPatchFieldMapper.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
const labelList & faceMap() const
Return face map.
static const zero Zero
Definition: zero.H:97
This boundary condition provides a waveVelocity condition. This sets the velocity to that specified b...
const vectorField & cellCentres() const
static const waveSuperposition & New(const objectRegistry &db)
Return a reference to the wave model on the given database,.
virtual label size() const
Return size.
Definition: fvPatch.H:155
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
volScalarField scalarField(fieldObject, mesh)
dimensionedScalar pos0(const dimensionedScalar &ds)
labelList f(nPoints)
Post-processing mesh subset tool. Given the original mesh and the list of selected cells...
Definition: fvMeshSubset.H:73
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
waveAlphaFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const fvMesh & subMesh() const
Return reference to subset mesh.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
ListType subset(const UList< T > &select, const T &value, const ListType &)
Extract elements of List when select is a certain value.
tmp< scalarField > alphan() const
Return the current modelled phase fraction field in the.
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.
tmp< vectorField > U() const
Return the current modelled velocity field on the patch faces.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540