GeometricBoundaryField.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 
26 #include "emptyPolyPatch.H"
27 #include "commSchedule.H"
28 #include "globalMeshData.H"
29 #include "cyclicPolyPatch.H"
30 
31 template<class Type, template<class> class PatchField, class GeoMesh>
34 (
36  const dictionary& dict
37 )
38 {
39  // Clear the boundary field if already initialised
40  this->clear();
41 
42  this->setSize(bmesh_.size());
43 
44  if (debug)
45  {
47  }
48 
49 
50  label nUnset = this->size();
51 
52  // 1. Handle explicit patch names. Note that there can be only one explicit
53  // patch name since is key of dictionary.
54  forAllConstIter(dictionary, dict, iter)
55  {
56  if (iter().isDict() && !iter().keyword().isPattern())
57  {
58  label patchi = bmesh_.findPatchID(iter().keyword());
59 
60  if (patchi != -1)
61  {
62  this->set
63  (
64  patchi,
66  (
67  bmesh_[patchi],
68  field,
69  iter().dict()
70  )
71  );
72  nUnset--;
73  }
74  }
75  }
76 
77  if (nUnset == 0)
78  {
79  return;
80  }
81 
82 
83  // 2. Patch-groups. (using non-wild card entries of dictionaries)
84  // (patchnames already matched above)
85  // Note: in reverse order of entries in the dictionary (last
86  // patchGroups wins). This is so is consistent with dictionary wildcard
87  // behaviour
88  if (dict.size())
89  {
90  for
91  (
93  iter != dict.rend();
94  ++iter
95  )
96  {
97  const entry& e = iter();
98 
99  if (e.isDict() && !e.keyword().isPattern())
100  {
101  const labelList patchIDs = bmesh_.findIndices
102  (
103  e.keyword(),
104  true // use patchGroups
105  );
106 
107  forAll(patchIDs, i)
108  {
109  label patchi = patchIDs[i];
110 
111  if (!this->set(patchi))
112  {
113  this->set
114  (
115  patchi,
117  (
118  bmesh_[patchi],
119  field,
120  e.dict()
121  )
122  );
123  }
124  }
125  }
126  }
127  }
128 
129 
130  // 3. Wildcard patch overrides
131  forAll(bmesh_, patchi)
132  {
133  if (!this->set(patchi))
134  {
135  if (bmesh_[patchi].type() == emptyPolyPatch::typeName)
136  {
137  this->set
138  (
139  patchi,
141  (
142  emptyPolyPatch::typeName,
143  bmesh_[patchi],
144  field
145  )
146  );
147  }
148  else
149  {
150  bool found = dict.found(bmesh_[patchi].name());
151 
152  if (found)
153  {
154  this->set
155  (
156  patchi,
158  (
159  bmesh_[patchi],
160  field,
161  dict.subDict(bmesh_[patchi].name())
162  )
163  );
164  }
165  }
166  }
167  }
168 
169 
170  // Check for any unset patches
171  forAll(bmesh_, patchi)
172  {
173  if (!this->set(patchi))
174  {
175  if (bmesh_[patchi].type() == cyclicPolyPatch::typeName)
176  {
178  (
179  dict
180  ) << "Cannot find patchField entry for cyclic "
181  << bmesh_[patchi].name() << endl
182  << "Is your field uptodate with split cyclics?" << endl
183  << "Run foamUpgradeCyclics to convert mesh and fields"
184  << " to split cyclics." << exit(FatalIOError);
185  }
186  else
187  {
189  (
190  dict
191  ) << "Cannot find patchField entry for "
192  << bmesh_[patchi].name() << exit(FatalIOError);
193  }
194  }
195  }
196 }
197 
198 
199 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
200 
201 template<class Type, template<class> class PatchField, class GeoMesh>
203 Boundary
204 (
205  const BoundaryMesh& bmesh
206 )
207 :
208  FieldField<PatchField, Type>(bmesh.size()),
209  bmesh_(bmesh)
210 {}
211 
212 
213 template<class Type, template<class> class PatchField, class GeoMesh>
215 Boundary
216 (
217  const BoundaryMesh& bmesh,
218  const DimensionedField<Type, GeoMesh>& field,
219  const word& patchFieldType
220 )
221 :
222  FieldField<PatchField, Type>(bmesh.size()),
223  bmesh_(bmesh)
224 {
225  if (debug)
226  {
227  InfoInFunction << endl;
228  }
229 
230  forAll(bmesh_, patchi)
231  {
232  this->set
233  (
234  patchi,
236  (
237  patchFieldType,
238  bmesh_[patchi],
239  field
240  )
241  );
242  }
243 }
244 
245 
246 template<class Type, template<class> class PatchField, class GeoMesh>
248 Boundary
249 (
250  const BoundaryMesh& bmesh,
251  const DimensionedField<Type, GeoMesh>& field,
252  const wordList& patchFieldTypes,
253  const wordList& constraintTypes
254 )
255 :
256  FieldField<PatchField, Type>(bmesh.size()),
257  bmesh_(bmesh)
258 {
259  if (debug)
260  {
261  InfoInFunction << endl;
262  }
263 
264  if
265  (
266  patchFieldTypes.size() != this->size()
267  || (constraintTypes.size() && (constraintTypes.size() != this->size()))
268  )
269  {
271  << "Incorrect number of patch type specifications given" << nl
272  << " Number of patches in mesh = " << bmesh.size()
273  << " number of patch type specifications = "
274  << patchFieldTypes.size()
275  << abort(FatalError);
276  }
277 
278  if (constraintTypes.size())
279  {
280  forAll(bmesh_, patchi)
281  {
282  this->set
283  (
284  patchi,
286  (
287  patchFieldTypes[patchi],
288  constraintTypes[patchi],
289  bmesh_[patchi],
290  field
291  )
292  );
293  }
294  }
295  else
296  {
297  forAll(bmesh_, patchi)
298  {
299  this->set
300  (
301  patchi,
303  (
304  patchFieldTypes[patchi],
305  bmesh_[patchi],
306  field
307  )
308  );
309  }
310  }
311 }
312 
313 
314 template<class Type, template<class> class PatchField, class GeoMesh>
316 Boundary
317 (
318  const BoundaryMesh& bmesh,
319  const DimensionedField<Type, GeoMesh>& field,
320  const PtrList<PatchField<Type>>& ptfl
321 )
322 :
323  FieldField<PatchField, Type>(bmesh.size()),
324  bmesh_(bmesh)
325 {
326  if (debug)
327  {
328  InfoInFunction << endl;
329  }
330 
331  forAll(bmesh_, patchi)
332  {
333  this->set(patchi, ptfl[patchi].clone(field));
334  }
335 }
336 
337 
338 template<class Type, template<class> class PatchField, class GeoMesh>
340 Boundary
341 (
342  const DimensionedField<Type, GeoMesh>& field,
344  Boundary& btf
345 )
346 :
348  bmesh_(btf.bmesh_)
349 {
350  if (debug)
351  {
352  InfoInFunction << endl;
353  }
354 
355  forAll(bmesh_, patchi)
356  {
357  this->set(patchi, btf[patchi].clone(field));
358  }
359 }
360 
361 
362 template<class Type, template<class> class PatchField, class GeoMesh>
364 Boundary
365 (
367  Boundary& btf
368 )
369 :
371  bmesh_(btf.bmesh_)
372 {
373  if (debug)
374  {
375  InfoInFunction << endl;
376  }
377 }
378 
379 
380 template<class Type, template<class> class PatchField, class GeoMesh>
382 Boundary
383 (
384  const BoundaryMesh& bmesh,
385  const DimensionedField<Type, GeoMesh>& field,
386  const dictionary& dict
387 )
388 :
389  FieldField<PatchField, Type>(bmesh.size()),
390  bmesh_(bmesh)
391 {
392  readField(field, dict);
393 }
394 
395 
396 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
397 
398 template<class Type, template<class> class PatchField, class GeoMesh>
401 {
402  if (debug)
403  {
404  InfoInFunction << endl;
405  }
406 
407  forAll(*this, patchi)
408  {
409  this->operator[](patchi).updateCoeffs();
410  }
411 }
412 
413 
414 template<class Type, template<class> class PatchField, class GeoMesh>
417 {
418  if (debug)
419  {
420  InfoInFunction << endl;
421  }
422 
423  if
424  (
427  )
428  {
429  label nReq = Pstream::nRequests();
430 
431  forAll(*this, patchi)
432  {
433  this->operator[](patchi).initEvaluate(Pstream::defaultCommsType);
434  }
435 
436  // Block for any outstanding requests
437  if
438  (
441  )
442  {
443  Pstream::waitRequests(nReq);
444  }
445 
446  forAll(*this, patchi)
447  {
449  }
450  }
452  {
453  const lduSchedule& patchSchedule =
454  bmesh_.mesh().globalData().patchSchedule();
455 
456  forAll(patchSchedule, patchEvali)
457  {
458  if (patchSchedule[patchEvali].init)
459  {
460  this->operator[](patchSchedule[patchEvali].patch)
461  .initEvaluate(Pstream::scheduled);
462  }
463  else
464  {
465  this->operator[](patchSchedule[patchEvali].patch)
466  .evaluate(Pstream::scheduled);
467  }
468  }
469  }
470  else
471  {
473  << "Unsuported communications type "
475  << exit(FatalError);
476  }
477 }
478 
479 
480 template<class Type, template<class> class PatchField, class GeoMesh>
483 types() const
484 {
485  const FieldField<PatchField, Type>& pff = *this;
486 
487  wordList Types(pff.size());
488 
489  forAll(pff, patchi)
490  {
491  Types[patchi] = pff[patchi].type();
492  }
493 
494  return Types;
495 }
496 
497 
498 template<class Type, template<class> class PatchField, class GeoMesh>
502 {
504  BoundaryInternalField(*this);
505 
506  forAll(BoundaryInternalField, patchi)
507  {
508  BoundaryInternalField[patchi] ==
509  this->operator[](patchi).patchInternalField();
510  }
511 
512  return BoundaryInternalField;
513 }
514 
515 
516 template<class Type, template<class> class PatchField, class GeoMesh>
519 interfaces() const
520 {
522 
523  forAll(interfaces, patchi)
524  {
525  if (isA<LduInterfaceField<Type>>(this->operator[](patchi)))
526  {
527  interfaces.set
528  (
529  patchi,
531  (
532  this->operator[](patchi)
533  )
534  );
535  }
536  }
537 
538  return interfaces;
539 }
540 
541 
542 template<class Type, template<class> class PatchField, class GeoMesh>
546 {
548 
549  forAll(interfaces, patchi)
550  {
551  if (isA<lduInterfaceField>(this->operator[](patchi)))
552  {
553  interfaces.set
554  (
555  patchi,
556  &refCast<const lduInterfaceField>
557  (
558  this->operator[](patchi)
559  )
560  );
561  }
562  }
563 
564  return interfaces;
565 }
566 
567 
568 template<class Type, template<class> class PatchField, class GeoMesh>
570 writeEntry(const word& keyword, Ostream& os) const
571 {
572  os << keyword << nl << token::BEGIN_BLOCK << incrIndent << nl;
573 
574  forAll(*this, patchi)
575  {
576  os << indent << this->operator[](patchi).patch().name() << nl
577  << indent << token::BEGIN_BLOCK << nl
578  << incrIndent << this->operator[](patchi) << decrIndent
579  << indent << token::END_BLOCK << endl;
580  }
581 
582  os << decrIndent << token::END_BLOCK << endl;
583 
584  // Check state of IOstream
585  os.check
586  (
587  "GeometricField<Type, PatchField, GeoMesh>::Boundary::"
588  "writeEntry(const word& keyword, Ostream& os) const"
589  );
590 }
591 
592 
593 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
594 
595 template<class Type, template<class> class PatchField, class GeoMesh>
597 operator=
598 (
600  Boundary& bf
601 )
602 {
604 }
605 
606 
607 template<class Type, template<class> class PatchField, class GeoMesh>
609 operator=
610 (
611  const FieldField<PatchField, Type>& ptff
612 )
613 {
615 }
616 
617 
618 template<class Type, template<class> class PatchField, class GeoMesh>
620 operator=
621 (
622  const Type& t
623 )
624 {
626 }
627 
628 
629 template<class Type, template<class> class PatchField, class GeoMesh>
631 operator==
632 (
634  Boundary& bf
635 )
636 {
637  forAll((*this), patchi)
638  {
639  this->operator[](patchi) == bf[patchi];
640  }
641 }
642 
643 
644 template<class Type, template<class> class PatchField, class GeoMesh>
646 operator==
647 (
648  const FieldField<PatchField, Type>& ptff
649 )
650 {
651  forAll((*this), patchi)
652  {
653  this->operator[](patchi) == ptff[patchi];
654  }
655 }
656 
657 
658 template<class Type, template<class> class PatchField, class GeoMesh>
660 operator==
661 (
662  const Type& t
663 )
664 {
665  forAll((*this), patchi)
666  {
667  this->operator[](patchi) == t;
668  }
669 }
670 
671 
672 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
673 
674 template<class Type, template<class> class PatchField, class GeoMesh>
675 Foam::Ostream& Foam::operator<<
676 (
677  Ostream& os,
679  Boundary& bf
680 )
681 {
682  os << static_cast<const FieldField<PatchField, Type>&>(bf);
683  return os;
684 }
685 
686 
687 // ************************************************************************* //
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
Boundary(const BoundaryMesh &)
Construct from a BoundaryMesh.
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 & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:223
GeoMesh::BoundaryMesh BoundaryMesh
Type of boundary mesh on which this.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
tUEqn clear()
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:106
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:633
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
void evaluate()
Evaluate boundary conditions.
An abstract base class for implicitly-coupled interface fields e.g. processor and cyclic patch fields...
virtual const dictionary & dict() const =0
Return dictionary if this entry is a dictionary.
static label nRequests()
Get number of outstanding requests.
Definition: UPstream.C:107
Generic GeometricField class.
wordList types() const
Return a list of the patch types.
Boundary boundaryInternalField() const
Return BoundaryField of the cell values neighbouring.
label size() const
Return number of elements in list.
Definition: DLListBaseI.H:77
lduInterfaceFieldPtrsList scalarInterfaces() const
Return a list of pointers for each patch field with only those.
void writeEntry(const word &keyword, Ostream &os) const
Write boundary field as dictionary entry.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
void readField(const Internal &field, const dictionary &dict)
Read the boundary field.
void operator=(const FieldField< Field, Type > &)
Definition: FieldField.C:284
void updateCoeffs()
Update the boundary condition coefficients.
bool set(const label) const
Is element set.
points setSize(newPointi)
const dimensionedScalar e
Elementary charge.
Definition: doubleFloat.H:78
A class for handling words, derived from string.
Definition: word.H:59
const keyType & keyword() const
Return keyword.
Definition: entry.H:120
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
const const_reverse_iterator & rend() const
Definition: UILList.H:347
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
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
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:237
static const NamedEnum< commsTypes, 3 > commsTypeNames
Definition: UPstream.H:71
static void waitRequests(const label start=0)
Wait until all requests (from start onwards) have finished.
Definition: UPstream.C:117
Intrusive doubly-linked list.
Definition: IDLList.H:47
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const T & operator[](const label) const
Return element const reference.
Definition: UPtrListI.H:103
static commsTypes defaultCommsType
Default commsType.
Definition: UPstream.H:268
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:393
label patchi
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
LduInterfaceFieldPtrsList< Type > interfaces() const
Return a list of pointers for each patch field with only those.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const_reverse_iterator rbegin() const
Definition: UILList.H:342
friend Ostream & operator(Ostream &, const FieldField< PatchField, Type > &)
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:230
tmp< DimensionedField< Type, GeoMesh > > clone() const
Clone.
bool found
virtual bool isDict() const
Return true if this entry is a dictionary.
Definition: entry.H:153
bool isPattern() const
Should be treated as a match rather than a literal string.
Definition: keyTypeI.H:76
A keyword and a list of tokens is an &#39;entry&#39;.
Definition: entry.H:65
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
IOerror FatalIOError
#define InfoInFunction
Report an information message using Foam::Info.