activeBaffleVelocityFvPatchVectorField.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) 2011-2019 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 "volFields.H"
29 #include "surfaceFields.H"
30 #include "cyclicFvPatch.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
36 (
37  const fvPatch& p,
39 )
40 :
41  fixedValueFvPatchVectorField(p, iF),
42  pName_("p"),
43  cyclicPatchName_(),
44  cyclicPatchLabel_(-1),
45  orientation_(1),
46  initWallSf_(0),
47  initCyclicSf_(0),
48  nbrCyclicSf_(0),
49  openFraction_(0),
50  openingTime_(0),
51  maxOpenFractionDelta_(0),
52  curTimeIndex_(-1)
53 {}
54 
55 
58 (
60  const fvPatch& p,
62  const fvPatchFieldMapper& mapper
63 )
64 :
65  fixedValueFvPatchVectorField(ptf, p, iF, mapper),
66  pName_(ptf.pName_),
67  cyclicPatchName_(ptf.cyclicPatchName_),
68  cyclicPatchLabel_(ptf.cyclicPatchLabel_),
69  orientation_(ptf.orientation_),
70  initWallSf_(ptf.initWallSf_),
71  initCyclicSf_(ptf.initCyclicSf_),
72  nbrCyclicSf_(ptf.nbrCyclicSf_),
73  openFraction_(ptf.openFraction_),
74  openingTime_(ptf.openingTime_),
75  maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
76  curTimeIndex_(-1)
77 {}
78 
79 
82 (
83  const fvPatch& p,
85  const dictionary& dict
86 )
87 :
88  fixedValueFvPatchVectorField(p, iF, dict, false),
89  pName_(dict.lookupOrDefault<word>("p", "p")),
90  cyclicPatchName_(dict.lookup("cyclicPatch")),
91  cyclicPatchLabel_(p.patch().boundaryMesh().findPatchID(cyclicPatchName_)),
92  orientation_(dict.lookup<label>("orientation")),
93  initWallSf_(p.Sf()),
94  initCyclicSf_(p.boundaryMesh()[cyclicPatchLabel_].Sf()),
95  nbrCyclicSf_
96  (
97  refCast<const cyclicFvPatch>
98  (
99  p.boundaryMesh()[cyclicPatchLabel_]
100  ).neighbFvPatch().Sf()
101  ),
102  openFraction_(dict.lookup<scalar>("openFraction")),
103  openingTime_(dict.lookup<scalar>("openingTime")),
104  maxOpenFractionDelta_(dict.lookup<scalar>("maxOpenFractionDelta")),
105  curTimeIndex_(-1)
106 {
107  fvPatchVectorField::operator=(Zero);
108 }
109 
110 
113 (
115 )
116 :
117  fixedValueFvPatchVectorField(ptf),
118  pName_(ptf.pName_),
119  cyclicPatchName_(ptf.cyclicPatchName_),
120  cyclicPatchLabel_(ptf.cyclicPatchLabel_),
121  orientation_(ptf.orientation_),
122  initWallSf_(ptf.initWallSf_),
123  initCyclicSf_(ptf.initCyclicSf_),
124  nbrCyclicSf_(ptf.nbrCyclicSf_),
125  openFraction_(ptf.openFraction_),
126  openingTime_(ptf.openingTime_),
127  maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
128  curTimeIndex_(-1)
129 {}
130 
131 
134 (
137 )
138 :
139  fixedValueFvPatchVectorField(ptf, iF),
140  pName_(ptf.pName_),
141  cyclicPatchName_(ptf.cyclicPatchName_),
142  cyclicPatchLabel_(ptf.cyclicPatchLabel_),
143  orientation_(ptf.orientation_),
144  initWallSf_(ptf.initWallSf_),
145  initCyclicSf_(ptf.initCyclicSf_),
146  nbrCyclicSf_(ptf.nbrCyclicSf_),
147  openFraction_(ptf.openFraction_),
148  openingTime_(ptf.openingTime_),
149  maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
150  curTimeIndex_(-1)
151 {}
152 
153 
154 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
155 
157 (
158  const fvPatchFieldMapper& m
159 )
160 {
161  fixedValueFvPatchVectorField::autoMap(m);
162 
163  //- Note: cannot map field from cyclic patch anyway so just recalculate
164  // Areas should be consistent when doing autoMap except in case of
165  // topo changes.
166  //- Note: we don't want to use Sf here since triggers rebuilding of
167  // fvMesh::S() which will give problems when mapped (since already
168  // on new mesh)
169  const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
170  initWallSf_ = patch().patchSlice(areas);
171  initCyclicSf_ = patch().boundaryMesh()
172  [
173  cyclicPatchLabel_
174  ].patchSlice(areas);
175  nbrCyclicSf_ = refCast<const cyclicFvPatch>
176  (
177  patch().boundaryMesh()
178  [
179  cyclicPatchLabel_
180  ]
181  ).neighbFvPatch().patch().patchSlice(areas);
182 }
183 
184 
186 (
187  const fvPatchVectorField& ptf,
188  const labelList& addr
189 )
190 {
191  fixedValueFvPatchVectorField::rmap(ptf, addr);
192 
193  // See autoMap.
194  const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
195  initWallSf_ = patch().patchSlice(areas);
196  initCyclicSf_ = patch().boundaryMesh()
197  [
198  cyclicPatchLabel_
199  ].patchSlice(areas);
200  nbrCyclicSf_ = refCast<const cyclicFvPatch>
201  (
202  patch().boundaryMesh()
203  [
204  cyclicPatchLabel_
205  ]
206  ).neighbFvPatch().patch().patchSlice(areas);
207 }
208 
209 
211 {
212  if (updated())
213  {
214  return;
215  }
216 
217  // Execute the change to the openFraction only once per time-step
218  if (curTimeIndex_ != this->db().time().timeIndex())
219  {
220  const volScalarField& p = db().lookupObject<volScalarField>
221  (
222  pName_
223  );
224 
225  const fvPatch& cyclicPatch = patch().boundaryMesh()[cyclicPatchLabel_];
226  const labelList& cyclicFaceCells = cyclicPatch.patch().faceCells();
227  const fvPatch& nbrPatch = refCast<const cyclicFvPatch>
228  (
229  cyclicPatch
230  ).neighbFvPatch();
231  const labelList& nbrFaceCells = nbrPatch.patch().faceCells();
232 
233  scalar forceDiff = 0;
234 
235  // Add this side
236  forAll(cyclicFaceCells, facei)
237  {
238  forceDiff += p[cyclicFaceCells[facei]]*mag(initCyclicSf_[facei]);
239  }
240 
241  // Remove other side
242  forAll(nbrFaceCells, facei)
243  {
244  forceDiff -= p[nbrFaceCells[facei]]*mag(nbrCyclicSf_[facei]);
245  }
246 
247  openFraction_ =
248  max
249  (
250  min
251  (
252  openFraction_
253  + min
254  (
255  this->db().time().deltaTValue()/openingTime_,
256  maxOpenFractionDelta_
257  )
258  *(orientation_*sign(forceDiff)),
259  1 - 1e-6
260  ),
261  1e-6
262  );
263 
264  Info<< "openFraction = " << openFraction_ << endl;
265 
266  vectorField::subField Sfw = this->patch().patch().faceAreas();
267  const vectorField newSfw((1 - openFraction_)*initWallSf_);
268  forAll(Sfw, facei)
269  {
270  Sfw[facei] = newSfw[facei];
271  }
272  const_cast<scalarField&>(patch().magSf()) = mag(patch().Sf());
273 
274  // Update owner side of cyclic
275  const_cast<vectorField&>(cyclicPatch.Sf()) =
276  openFraction_*initCyclicSf_;
277  const_cast<scalarField&>(cyclicPatch.magSf()) =
278  mag(cyclicPatch.Sf());
279  // Update neighbour side of cyclic
280  const_cast<vectorField&>(nbrPatch.Sf()) =
281  openFraction_*nbrCyclicSf_;
282  const_cast<scalarField&>(nbrPatch.magSf()) =
283  mag(nbrPatch.Sf());
284 
285  curTimeIndex_ = this->db().time().timeIndex();
286  }
287 
288  fixedValueFvPatchVectorField::updateCoeffs();
289 }
290 
291 
293 {
295  writeEntryIfDifferent<word>(os, "p", "p", pName_);
296  writeEntry(os, "cyclicPatch", cyclicPatchName_);
297  writeEntry(os, "orientation", orientation_);
298  writeEntry(os, "openingTime", openingTime_);
299  writeEntry(os, "maxOpenFractionDelta", maxOpenFractionDelta_);
300  writeEntry(os, "openFraction", openFraction_);
301  writeEntryIfDifferent<word>(os, "p", "p", pName_);
302  writeEntry(os, "value", *this);
303 }
304 
305 
306 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
307 
308 namespace Foam
309 {
311  (
314  );
315 }
316 
317 
318 // ************************************************************************* //
Foam::surfaceFields.
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:179
dimensionedScalar sign(const dimensionedScalar &ds)
#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
const scalarField & magSf() const
Return face area magnitudes.
Definition: fvPatch.C:136
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:278
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
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
label findPatchID(const word &patchName) const
Find patch index given a name.
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:280
Pre-declare related SubField type.
Definition: Field.H:60
Macros for easy insertion into run-time selection tables.
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:321
A class for handling words, derived from string.
Definition: word.H:59
activeBaffleVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:137
virtual void rmap(const fvPatchVectorField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Foam::fvPatchFieldMapper.
static const zero Zero
Definition: zero.H:97
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const vectorField & Sf() const
Return face area vectors.
Definition: fvPatch.C:130
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
label timeIndex
Definition: getTimeIndex.H:4
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
This velocity boundary condition simulates the opening of a baffle due to local flow conditions...
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812