activePressureForceBaffleVelocityFvPatchVectorField.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-2016 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  minThresholdValue_(0),
54  fBased_(1),
55  baffleActivated_(0)
56 {}
57 
58 
61 (
63  const fvPatch& p,
65  const fvPatchFieldMapper& mapper
66 )
67 :
68  fixedValueFvPatchVectorField(ptf, p, iF, mapper),
69  pName_(ptf.pName_),
70  cyclicPatchName_(ptf.cyclicPatchName_),
71  cyclicPatchLabel_(ptf.cyclicPatchLabel_),
72  orientation_(ptf.orientation_),
73  initWallSf_(ptf.initWallSf_),
74  initCyclicSf_(ptf.initCyclicSf_),
75  nbrCyclicSf_(ptf.nbrCyclicSf_),
76  openFraction_(ptf.openFraction_),
77  openingTime_(ptf.openingTime_),
78  maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
79  curTimeIndex_(-1),
80  minThresholdValue_(ptf.minThresholdValue_),
81  fBased_(ptf.fBased_),
82  baffleActivated_(ptf.baffleActivated_)
83 {}
84 
85 
88 (
89  const fvPatch& p,
91  const dictionary& dict
92 )
93 :
94  fixedValueFvPatchVectorField(p, iF),
95  pName_(dict.lookupOrDefault<word>("p", "p")),
96  cyclicPatchName_(dict.lookup("cyclicPatch")),
97  cyclicPatchLabel_(p.patch().boundaryMesh().findPatchID(cyclicPatchName_)),
98  orientation_(readLabel(dict.lookup("orientation"))),
99  initWallSf_(0),
100  initCyclicSf_(0),
101  nbrCyclicSf_(0),
102  openFraction_(readScalar(dict.lookup("openFraction"))),
103  openingTime_(readScalar(dict.lookup("openingTime"))),
104  maxOpenFractionDelta_(readScalar(dict.lookup("maxOpenFractionDelta"))),
105  curTimeIndex_(-1),
106  minThresholdValue_(readScalar(dict.lookup("minThresholdValue"))),
107  fBased_(readBool(dict.lookup("forceBased"))),
108  baffleActivated_(0)
109 {
110  fvPatchVectorField::operator=(Zero);
111 
112  if (p.size() > 0)
113  {
114  initWallSf_ = p.Sf();
115  initCyclicSf_ = p.boundaryMesh()[cyclicPatchLabel_].Sf();
116  nbrCyclicSf_ = refCast<const cyclicFvPatch>
117  (
118  p.boundaryMesh()[cyclicPatchLabel_]
119  ).neighbFvPatch().Sf();
120  }
121 
122  if (dict.found("p"))
123  {
124  dict.lookup("p") >> pName_;
125  }
126 }
127 
128 
131 (
133 )
134 :
135  fixedValueFvPatchVectorField(ptf),
136  pName_(ptf.pName_),
137  cyclicPatchName_(ptf.cyclicPatchName_),
138  cyclicPatchLabel_(ptf.cyclicPatchLabel_),
139  orientation_(ptf.orientation_),
140  initWallSf_(ptf.initWallSf_),
141  initCyclicSf_(ptf.initCyclicSf_),
142  nbrCyclicSf_(ptf.nbrCyclicSf_),
143  openFraction_(ptf.openFraction_),
144  openingTime_(ptf.openingTime_),
145  maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
146  curTimeIndex_(-1),
147  minThresholdValue_(ptf.minThresholdValue_),
148  fBased_(ptf.fBased_),
149  baffleActivated_(ptf.baffleActivated_)
150 {}
151 
152 
155 (
158 )
159 :
160  fixedValueFvPatchVectorField(ptf, iF),
161  pName_(ptf.pName_),
162  cyclicPatchName_(ptf.cyclicPatchName_),
163  cyclicPatchLabel_(ptf.cyclicPatchLabel_),
164  orientation_(ptf.orientation_),
165  initWallSf_(ptf.initWallSf_),
166  initCyclicSf_(ptf.initCyclicSf_),
167  nbrCyclicSf_(ptf.nbrCyclicSf_),
168  openFraction_(ptf.openFraction_),
169  openingTime_(ptf.openingTime_),
170  maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
171  curTimeIndex_(-1),
172  minThresholdValue_(ptf.minThresholdValue_),
173  fBased_(ptf.fBased_),
174  baffleActivated_(ptf.baffleActivated_)
175 {}
176 
177 
178 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
179 
181 (
182  const fvPatchFieldMapper& m
183 )
184 {
185  fixedValueFvPatchVectorField::autoMap(m);
186 
187  //- Note: cannot map field from cyclic patch anyway so just recalculate
188  // Areas should be consistent when doing autoMap except in case of
189  // topo changes.
190  //- Note: we don't want to use Sf here since triggers rebuilding of
191  // fvMesh::S() which will give problems when mapped (since already
192  // on new mesh)
193  forAll(patch().boundaryMesh().mesh().faceAreas(), i)
194  {
195  if (mag(patch().boundaryMesh().mesh().faceAreas()[i]) == 0)
196  {
197  Info << "faceArea[active] "<< i << endl;
198  }
199  }
200  if (patch().size() > 0)
201  {
202  const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
203  initWallSf_ = patch().patchSlice(areas);
204  initCyclicSf_ = patch().boundaryMesh()
205  [
206  cyclicPatchLabel_
207  ].patchSlice(areas);
208  nbrCyclicSf_ = refCast<const cyclicFvPatch>
209  (
210  patch().boundaryMesh()
211  [
212  cyclicPatchLabel_
213  ]
214  ).neighbFvPatch().patch().patchSlice(areas);
215  }
216 }
217 
219 (
220  const fvPatchVectorField& ptf,
221  const labelList& addr
222 )
223 {
224  fixedValueFvPatchVectorField::rmap(ptf, addr);
225 
226  // See autoMap.
227  const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
228  initWallSf_ = patch().patchSlice(areas);
229  initCyclicSf_ = patch().boundaryMesh()
230  [
231  cyclicPatchLabel_
232  ].patchSlice(areas);
233  nbrCyclicSf_ = refCast<const cyclicFvPatch>
234  (
235  patch().boundaryMesh()
236  [
237  cyclicPatchLabel_
238  ]
239  ).neighbFvPatch().patch().patchSlice(areas);
240 }
241 
242 
244 {
245  if (updated())
246  {
247  return;
248  }
249  // Execute the change to the openFraction only once per time-step
250  if (curTimeIndex_ != this->db().time().timeIndex())
251  {
252  const volScalarField& p = db().lookupObject<volScalarField>
253  (
254  pName_
255  );
256 
257  const fvPatch& cyclicPatch = patch().boundaryMesh()[cyclicPatchLabel_];
258  const labelList& cyclicFaceCells = cyclicPatch.patch().faceCells();
259  const fvPatch& nbrPatch = refCast<const cyclicFvPatch>
260  (
261  cyclicPatch
262  ).neighbFvPatch();
263 
264  const labelList& nbrFaceCells = nbrPatch.patch().faceCells();
265 
266  scalar valueDiff = 0;
267 
268  if (fBased_)
269  {
270  // Add this side
271  forAll(cyclicFaceCells, facei)
272  {
273  valueDiff +=p[cyclicFaceCells[facei]]*mag(initCyclicSf_[facei]);
274  }
275 
276  // Remove other side
277  forAll(nbrFaceCells, facei)
278  {
279  valueDiff -=p[nbrFaceCells[facei]]*mag(initCyclicSf_[facei]);
280  }
281 
282  Info<< "Force difference = " << valueDiff << endl;
283  }
284  else //pressure based
285  {
286  forAll(cyclicFaceCells, facei)
287  {
288  valueDiff += p[cyclicFaceCells[facei]];
289  }
290 
291  forAll(nbrFaceCells, facei)
292  {
293  valueDiff -= p[nbrFaceCells[facei]];
294  }
295 
296  Info<< "Pressure difference = " << valueDiff << endl;
297  }
298 
299  if ((mag(valueDiff) > mag(minThresholdValue_)) || baffleActivated_)
300  {
301  openFraction_ =
302  max(
303  min(
304  openFraction_
305  + min
306  (
307  this->db().time().deltaT().value()/openingTime_,
308  maxOpenFractionDelta_
309  )*(orientation_),
310  1 - 1e-6
311  ),
312  1e-6
313  );
314 
315  baffleActivated_ = true;
316  }
317  else
318  {
319  openFraction_ = max(min(1 - 1e-6, openFraction_), 1e-6);
320  }
321 
322  Info<< "Open fraction = " << openFraction_ << endl;
323 
324  vectorField::subField Sfw = patch().patch().faceAreas();
325  vectorField newSfw((1 - openFraction_)*initWallSf_);
326  forAll(Sfw, facei)
327  {
328  Sfw[facei] = newSfw[facei];
329  }
330  const_cast<scalarField&>(patch().magSf()) = mag(patch().Sf());
331 
332  // Update owner side of cyclic
333  const_cast<vectorField&>(cyclicPatch.Sf()) =
334  openFraction_*initCyclicSf_;
335 
336  const_cast<scalarField&>(cyclicPatch.magSf()) =
337  mag(cyclicPatch.Sf());
338 
339  // Update neighbour side of cyclic
340  const_cast<vectorField&>(nbrPatch.Sf()) =
341  openFraction_*nbrCyclicSf_;
342 
343  const_cast<scalarField&>(nbrPatch.magSf()) =
344  mag(nbrPatch.Sf());
345 
346  curTimeIndex_ = this->db().time().timeIndex();
347  }
348 
349  fixedValueFvPatchVectorField::updateCoeffs();
350 }
351 
352 
354 write(Ostream& os) const
355 {
357  writeEntryIfDifferent<word>(os, "p", "p", pName_);
358  os.writeKeyword("cyclicPatch")
359  << cyclicPatchName_ << token::END_STATEMENT << nl;
360  os.writeKeyword("orientation")
361  << orientation_ << token::END_STATEMENT << nl;
362  os.writeKeyword("openingTime")
363  << openingTime_ << token::END_STATEMENT << nl;
364  os.writeKeyword("maxOpenFractionDelta")
365  << maxOpenFractionDelta_ << token::END_STATEMENT << nl;
366  os.writeKeyword("openFraction")
367  << openFraction_ << token::END_STATEMENT << nl;
368  os.writeKeyword("minThresholdValue")
369  << minThresholdValue_ << token::END_STATEMENT << nl;
370  os.writeKeyword("forceBased")
371  << fBased_ << token::END_STATEMENT << nl;
372  writeEntry("value", os);
373 }
374 
375 
376 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
377 
378 namespace Foam
379 {
381  (
384  );
385 }
386 
387 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
Foam::surfaceFields.
This boundary condition is applied to the flow velocity, to simulate the opening or closure of a baff...
const scalarField & magSf() const
Return face area magnitudes.
Definition: fvPatch.C:136
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const double e
Elementary charge.
Definition: doubleFloat.H:78
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.
const vectorField & Sf() const
Return face area vectors.
Definition: fvPatch.C:130
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
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:65
bool readBool(Istream &)
Definition: boolIO.C:60
Pre-declare related SubField type.
Definition: Field.H:61
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:143
Macros for easy insertion into run-time selection tables.
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
activePressureForceBaffleVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
Foam::fvPatchFieldMapper.
virtual label size() const
Return size.
Definition: fvPatch.H:161
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
static const zero Zero
Definition: zero.H:91
label readLabel(Istream &is)
Definition: label.H:64
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:262
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:315
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 > &)
label findPatchID(const word &patchName) const
Find patch index given a name.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:278
virtual void rmap(const fvPatchVectorField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Namespace for OpenFOAM.
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:363
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:185