waveVelocityFvPatchVectorField.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 
29 #include "levelSet.H"
30 #include "volFields.H"
31 #include "fvMeshSubset.H"
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
36 (
37  const fvPatch& p,
39 )
40 :
41  directionMixedFvPatchVectorField(p, iF),
42  phiName_("phi"),
43  pName_("p"),
44  inletOutlet_(true),
45  faceCellSubset_(nullptr),
46  faceCellSubsetTimeIndex_(-1)
47 {
48  refValue() = Zero;
49  refGrad() = Zero;
50  valueFraction() = Zero;
51 }
52 
53 
55 (
56  const fvPatch& p,
58  const dictionary& dict
59 )
60 :
61  directionMixedFvPatchVectorField(p, iF),
62  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
63  pName_(dict.lookupOrDefault<word>("p", "p")),
64  inletOutlet_(dict.lookupOrDefault<Switch>("inletOutlet", true)),
65  faceCellSubset_(nullptr),
66  faceCellSubsetTimeIndex_(-1)
67 {
68  if (dict.found("value"))
69  {
70  fvPatchVectorField::operator=(vectorField("value", dict, p.size()));
71  }
72  else
73  {
74  fvPatchVectorField::operator=(patchInternalField());
75  }
76 
77  refValue() = *this;
78  refGrad() = Zero;
79  valueFraction() = Zero;
80 }
81 
82 
84 (
86  const fvPatch& p,
88  const fvPatchFieldMapper& mapper
89 )
90 :
91  directionMixedFvPatchVectorField(ptf, p, iF, mapper),
92  phiName_(ptf.phiName_),
93  pName_(ptf.pName_),
94  inletOutlet_(ptf.inletOutlet_),
95  faceCellSubset_(nullptr),
96  faceCellSubsetTimeIndex_(-1)
97 {}
98 
99 
101 (
103 )
104 :
105  directionMixedFvPatchVectorField(ptf),
106  phiName_(ptf.phiName_),
107  pName_(ptf.pName_),
108  inletOutlet_(ptf.inletOutlet_),
109  faceCellSubset_(nullptr),
110  faceCellSubsetTimeIndex_(-1)
111 {}
112 
113 
115 (
118 )
119 :
120  directionMixedFvPatchVectorField(ptf, iF),
121  phiName_(ptf.phiName_),
122  pName_(ptf.pName_),
123  inletOutlet_(ptf.inletOutlet_),
124  faceCellSubset_(nullptr),
125  faceCellSubsetTimeIndex_(-1)
126 {}
127 
128 
129 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
130 
131 const Foam::fvMeshSubset&
133 {
134  const fvMesh& mesh = patch().boundaryMesh().mesh();
135  const label timeIndex = mesh.time().timeIndex();
136 
137  if
138  (
139  !faceCellSubset_.valid()
140  || (mesh.changing() && faceCellSubsetTimeIndex_ != timeIndex)
141  )
142  {
143  faceCellSubset_.reset(new fvMeshSubset(mesh));
144  faceCellSubset_->setCellSubset(patch().faceCells());
145  faceCellSubsetTimeIndex_ = timeIndex;
146 
147  // Ask for the tetBasePtIs to trigger all processors to build them.
148  // Without this, processors that do not contain this patch will
149  // generate a comms mismatch.
150  faceCellSubset_->subMesh().tetBasePtIs();
151  }
152 
153  return faceCellSubset_();
154 }
155 
156 
158 {
159  const scalar t = db().time().timeOutputValue();
160  const waveSuperposition& waves = waveSuperposition::New(db());
161 
162  return
164  (
165  patch(),
166  waves.height(t, patch().Cf()),
167  waves.height(t, patch().patch().localPoints()),
168  waves.UGas(t, patch().Cf())(),
169  waves.UGas(t, patch().patch().localPoints())(),
170  waves.ULiquid(t, patch().Cf())(),
171  waves.ULiquid(t, patch().patch().localPoints())()
172  );
173 }
174 
175 
177 {
178  const scalar t = db().time().timeOutputValue();
179  const waveSuperposition& waves = waveSuperposition::New(db());
180 
182  const fvMesh& meshs = subset.subMesh();
183  const label patchis = findIndex(subset.patchMap(), patch().index());
184 
185  const vectorField Us
186  (
188  (
189  meshs,
190  waves.height(t, meshs.cellCentres())(),
191  waves.height(t, meshs.points())(),
192  waves.UGas(t, meshs.cellCentres())(),
193  waves.UGas(t, meshs.points())(),
194  waves.ULiquid(t, meshs.cellCentres())(),
195  waves.ULiquid(t, meshs.points())()
196  )
197  );
198 
199  tmp<vectorField> tResult(new vectorField(patch().size()));
200  vectorField& result = tResult.ref();
201 
202  if (patchis != -1)
203  {
204  forAll(meshs.boundary()[patchis], is)
205  {
206  const label fs = is + meshs.boundary()[patchis].patch().start();
207  const label cs = meshs.boundary()[patchis].faceCells()[is];
208  const label f = subset.faceMap()[fs];
209  const label i = patch().patch().whichFace(f);
210  result[i] = Us[cs];
211  }
212  }
213 
214  return tResult;
215 }
216 
217 
219 {
220  if (updated())
221  {
222  return;
223  }
224 
225  const fvPatchScalarField& pp =
226  patch().lookupPatchField<volScalarField, scalar>(pName_);
227 
228  if (isA<wavePressureFvPatchScalarField>(pp))
229  {
230  const vectorField U(this->U()), Un(this->Un());
231  const scalarField out(pos0(U & patch().Sf()));
232 
233  // Where inflow, set all velocity components to values specified by the
234  // wave model. Where outflow, set the tangential values and the normal
235  // gradient.
236  valueFraction() = symmTensor::I - out*sqr(patch().nf());
237  refValue() = U;
238  refGrad() = (U - Un)*patch().deltaCoeffs();
239  }
240  else
241  {
242  const vectorField U(this->U());
243 
244  if (inletOutlet_)
245  {
246  const scalarField& phip =
247  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
248  const scalarField out(pos0(phip));
249 
250  // Where inflow, fix all velocity components to values specified by
251  // the wave model.
252  refValue() = (1 - out)*U;
253  valueFraction() = (1 - out)*symmTensor::I;
254 
255  // Where outflow, set the normal component of the velocity to a
256  // value consistent with phi, but scale it to get the volumetric
257  // flow rate specified by the wave model. Tangential components are
258  // extrapolated.
259  const scalar QPhip = gSum(out*phip);
260  const scalar QWave = gSum(out*(U & patch().Sf()));
261  const vectorField nBySf(patch().Sf()/sqr(patch().magSf()));
262  if (QPhip > vSmall)
263  {
264  refValue() += out*(QWave/QPhip)*phip*nBySf;
265  }
266  else
267  {
268  refValue() += out*QWave*nBySf;
269  }
270  valueFraction() += out*sqr(patch().nf());
271  }
272  else
273  {
274  refValue() = U;
275  valueFraction() = symmTensor::I;
276  }
277  }
278 
279  directionMixedFvPatchVectorField::updateCoeffs();
280  directionMixedFvPatchVectorField::evaluate();
281 }
282 
283 
285 (
286  Ostream& os
287 ) const
288 {
290  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
291  writeEntryIfDifferent<word>(os, "p", "p", pName_);
292  writeEntryIfDifferent<Switch>(os, "inletOutlet", true, inletOutlet_);
293 }
294 
295 
296 // * * * * * * * * * * * * * * Build Macro Function * * * * * * * * * * * * //
297 
298 namespace Foam
299 {
301  (
304  );
305 }
306 
307 // ************************************************************************* //
bool changing() const
Is mesh changing (topology changing and/or moving)
Definition: polyMesh.H:540
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:438
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual tmp< vectorField > UGas(const scalar t, const vectorField &p) const
Get the gas velocity at a given time and global positions.
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
waveVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
const fvMeshSubset & faceCellSubset() const
Access the face-cell subset.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual void write(Ostream &) const
Write.
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
virtual tmp< vectorField > ULiquid(const scalar t, const vectorField &p) const
Get the liquid velocity at a given time and global positions.
volVectorField vectorField(fieldObject, mesh)
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
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:52
void write(Ostream &, const label, const dictionary &)
Write with dictionary lookup.
virtual tmp< scalarField > height(const scalar t, const vectorField &p) const
Get the height above the waves at a given time and global positions.
dynamicFvMesh & mesh
Type gSum(const FieldField< Field, Type > &f)
static const SymmTensor I
Definition: SymmTensor.H:72
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
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:53
dimensionedScalar pos0(const dimensionedScalar &ds)
tmp< Field< Type > > levelSetAverage(const fvMesh &mesh, const scalarField &levelC, const scalarField &levelP, const Field< Type > &positiveC, const Field< Type > &positiveP, const Field< Type > &negativeC, const Field< Type > &negativeP)
Calculate the average value of two fields, one on each side of a level set.
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,.
tmp< vectorField > Un() const
Return the current modelled velocity field in the neighbour cell.
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
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
label timeIndex
Definition: getTimeIndex.H:4
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.
Field< vector > vectorField
Specialisation of Field<T> for vector.
mesh Sf()
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