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