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-2013 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  {
46  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
47  "GeometricBoundaryField::readField"
48  "("
49  "const DimensionedField<Type, GeoMesh>&, "
50  "const dictionary&"
51  ")"
52  << endl;
53  }
54 
55 
56  label nUnset = this->size();
57 
58  // 1. Handle explicit patch names. Note that there can be only one explicit
59  // patch name since is key of dictionary.
60  forAllConstIter(dictionary, dict, iter)
61  {
62  if (iter().isDict() && !iter().keyword().isPattern())
63  {
64  label patchi = bmesh_.findPatchID(iter().keyword());
65 
66  if (patchi != -1)
67  {
68  this->set
69  (
70  patchi,
72  (
73  bmesh_[patchi],
74  field,
75  iter().dict()
76  )
77  );
78  nUnset--;
79  }
80  }
81  }
82 
83  if (nUnset == 0)
84  {
85  return;
86  }
87 
88 
89  // 2. Patch-groups. (using non-wild card entries of dictionaries)
90  // (patchnames already matched above)
91  // Note: in reverse order of entries in the dictionary (last
92  // patchGroups wins). This is so is consistent with dictionary wildcard
93  // behaviour
94  if (dict.size())
95  {
96  for
97  (
99  iter != dict.rend();
100  ++iter
101  )
102  {
103  const entry& e = iter();
104 
105  if (e.isDict() && !e.keyword().isPattern())
106  {
107  const labelList patchIDs = bmesh_.findIndices
108  (
109  e.keyword(),
110  true // use patchGroups
111  );
112 
113  forAll(patchIDs, i)
114  {
115  label patchi = patchIDs[i];
116 
117  if (!this->set(patchi))
118  {
119  this->set
120  (
121  patchi,
123  (
124  bmesh_[patchi],
125  field,
126  e.dict()
127  )
128  );
129  }
130  }
131  }
132  }
133  }
134 
135 
136  // 3. Wildcard patch overrides
137  forAll(bmesh_, patchi)
138  {
139  if (!this->set(patchi))
140  {
141  if (bmesh_[patchi].type() == emptyPolyPatch::typeName)
142  {
143  this->set
144  (
145  patchi,
147  (
148  emptyPolyPatch::typeName,
149  bmesh_[patchi],
150  field
151  )
152  );
153  }
154  else
155  {
156  bool found = dict.found(bmesh_[patchi].name());
157 
158  if (found)
159  {
160  this->set
161  (
162  patchi,
164  (
165  bmesh_[patchi],
166  field,
167  dict.subDict(bmesh_[patchi].name())
168  )
169  );
170  }
171  }
172  }
173  }
174 
175 
176  // Check for any unset patches
177  forAll(bmesh_, patchi)
178  {
179  if (!this->set(patchi))
180  {
181  if (bmesh_[patchi].type() == cyclicPolyPatch::typeName)
182  {
184  (
185  "GeometricField<Type, PatchField, GeoMesh>::"
186  "GeometricBoundaryField::readField"
187  "("
188  "const DimensionedField<Type, GeoMesh>&, "
189  "const dictionary&"
190  ")",
191  dict
192  ) << "Cannot find patchField entry for cyclic "
193  << bmesh_[patchi].name() << endl
194  << "Is your field uptodate with split cyclics?" << endl
195  << "Run foamUpgradeCyclics to convert mesh and fields"
196  << " to split cyclics." << exit(FatalIOError);
197  }
198  else
199  {
201  (
202  "GeometricField<Type, PatchField, GeoMesh>::"
203  "GeometricBoundaryField::readField"
204  "("
205  "const DimensionedField<Type, GeoMesh>&, "
206  "const dictionary&"
207  ")",
208  dict
209  ) << "Cannot find patchField entry for "
210  << bmesh_[patchi].name() << exit(FatalIOError);
211  }
212  }
213  }
214 }
215 
216 
217 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
218 
219 template<class Type, template<class> class PatchField, class GeoMesh>
222 (
223  const BoundaryMesh& bmesh
224 )
225 :
226  FieldField<PatchField, Type>(bmesh.size()),
227  bmesh_(bmesh)
228 {}
229 
230 
231 template<class Type, template<class> class PatchField, class GeoMesh>
234 (
235  const BoundaryMesh& bmesh,
236  const DimensionedField<Type, GeoMesh>& field,
237  const word& patchFieldType
238 )
239 :
240  FieldField<PatchField, Type>(bmesh.size()),
241  bmesh_(bmesh)
242 {
243  if (debug)
244  {
245  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
246  "GeometricBoundaryField::"
247  "GeometricBoundaryField(const BoundaryMesh&, "
248  "const DimensionedField<Type>&, const word&)"
249  << endl;
250  }
251 
252  forAll(bmesh_, patchi)
253  {
254  this->set
255  (
256  patchi,
258  (
259  patchFieldType,
260  bmesh_[patchi],
261  field
262  )
263  );
264  }
265 }
266 
267 
268 template<class Type, template<class> class PatchField, class GeoMesh>
271 (
272  const BoundaryMesh& bmesh,
273  const DimensionedField<Type, GeoMesh>& field,
274  const wordList& patchFieldTypes,
275  const wordList& constraintTypes
276 )
277 :
278  FieldField<PatchField, Type>(bmesh.size()),
279  bmesh_(bmesh)
280 {
281  if (debug)
282  {
283  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
284  "GeometricBoundaryField::"
285  "GeometricBoundaryField"
286  "("
287  "const BoundaryMesh&, "
288  "const DimensionedField<Type>&, "
289  "const wordList&, "
290  "const wordList&"
291  ")"
292  << endl;
293  }
294 
295  if
296  (
297  patchFieldTypes.size() != this->size()
298  || (constraintTypes.size() && (constraintTypes.size() != this->size()))
299  )
300  {
302  (
303  "GeometricField<Type, PatchField, GeoMesh>::"
304  "GeometricBoundaryField::"
305  "GeometricBoundaryField"
306  "("
307  "const BoundaryMesh&, "
308  "const DimensionedField<Type>&, "
309  "const wordList&, "
310  "const wordList&"
311  ")"
312  ) << "Incorrect number of patch type specifications given" << nl
313  << " Number of patches in mesh = " << bmesh.size()
314  << " number of patch type specifications = "
315  << patchFieldTypes.size()
316  << abort(FatalError);
317  }
318 
319  if (constraintTypes.size())
320  {
321  forAll(bmesh_, patchi)
322  {
323  this->set
324  (
325  patchi,
327  (
328  patchFieldTypes[patchi],
329  constraintTypes[patchi],
330  bmesh_[patchi],
331  field
332  )
333  );
334  }
335  }
336  else
337  {
338  forAll(bmesh_, patchi)
339  {
340  this->set
341  (
342  patchi,
344  (
345  patchFieldTypes[patchi],
346  bmesh_[patchi],
347  field
348  )
349  );
350  }
351  }
352 }
353 
354 
355 template<class Type, template<class> class PatchField, class GeoMesh>
358 (
359  const BoundaryMesh& bmesh,
360  const DimensionedField<Type, GeoMesh>& field,
361  const PtrList<PatchField<Type> >& ptfl
362 )
363 :
364  FieldField<PatchField, Type>(bmesh.size()),
365  bmesh_(bmesh)
366 {
367  if (debug)
368  {
369  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
370  "GeometricBoundaryField::"
371  "GeometricBoundaryField"
372  "("
373  "const BoundaryMesh&, "
374  "const DimensionedField<Type, GeoMesh>&, "
375  "const PtrLIst<PatchField<Type> >&"
376  ")"
377  << endl;
378  }
379 
380  forAll(bmesh_, patchi)
381  {
382  this->set(patchi, ptfl[patchi].clone(field));
383  }
384 }
385 
386 
387 template<class Type, template<class> class PatchField, class GeoMesh>
390 (
391  const DimensionedField<Type, GeoMesh>& field,
394 )
395 :
397  bmesh_(btf.bmesh_)
398 {
399  if (debug)
400  {
401  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
402  "GeometricBoundaryField::"
403  "GeometricBoundaryField"
404  "("
405  "const DimensionedField<Type, GeoMesh>&, "
406  "const typename GeometricField<Type, PatchField, GeoMesh>::"
407  "GeometricBoundaryField&"
408  ")"
409  << endl;
410  }
411 
412  forAll(bmesh_, patchi)
413  {
414  this->set(patchi, btf[patchi].clone(field));
415  }
416 }
417 
418 
419 // Construct as copy
420 // Dangerous because Field may be set to a field which gets deleted.
421 // Need new type of GeometricBoundaryField, one which IS part of a geometric
422 // field for which snGrad etc. may be called and a free standing
423 // GeometricBoundaryField for which such operations are unavailable.
424 template<class Type, template<class> class PatchField, class GeoMesh>
427 (
430 )
431 :
433  bmesh_(btf.bmesh_)
434 {
435  if (debug)
436  {
437  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
438  "GeometricBoundaryField::"
439  "GeometricBoundaryField"
440  "("
441  "const GeometricField<Type, PatchField, GeoMesh>::"
442  "GeometricBoundaryField&"
443  ")"
444  << endl;
445  }
446 }
447 
448 
449 template<class Type, template<class> class PatchField, class GeoMesh>
452 (
453  const BoundaryMesh& bmesh,
454  const DimensionedField<Type, GeoMesh>& field,
455  const dictionary& dict
456 )
457 :
458  FieldField<PatchField, Type>(bmesh.size()),
459  bmesh_(bmesh)
460 {
461  readField(field, dict);
462 }
463 
464 
465 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
466 
467 template<class Type, template<class> class PatchField, class GeoMesh>
470 {
471  if (debug)
472  {
473  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
474  "GeometricBoundaryField::"
475  "updateCoeffs()" << endl;
476  }
477 
478  forAll(*this, patchi)
479  {
480  this->operator[](patchi).updateCoeffs();
481  }
482 }
483 
484 
485 template<class Type, template<class> class PatchField, class GeoMesh>
488 {
489  if (debug)
490  {
491  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
492  "GeometricBoundaryField::"
493  "evaluate()" << endl;
494  }
495 
496  if
497  (
500  )
501  {
502  label nReq = Pstream::nRequests();
503 
504  forAll(*this, patchi)
505  {
506  this->operator[](patchi).initEvaluate(Pstream::defaultCommsType);
507  }
508 
509  // Block for any outstanding requests
510  if
511  (
514  )
515  {
516  Pstream::waitRequests(nReq);
517  }
518 
519  forAll(*this, patchi)
520  {
522  }
523  }
525  {
526  const lduSchedule& patchSchedule =
527  bmesh_.mesh().globalData().patchSchedule();
528 
529  forAll(patchSchedule, patchEvali)
530  {
531  if (patchSchedule[patchEvali].init)
532  {
533  this->operator[](patchSchedule[patchEvali].patch)
534  .initEvaluate(Pstream::scheduled);
535  }
536  else
537  {
538  this->operator[](patchSchedule[patchEvali].patch)
539  .evaluate(Pstream::scheduled);
540  }
541  }
542  }
543  else
544  {
545  FatalErrorIn("GeometricBoundaryField::evaluate()")
546  << "Unsuported communications type "
548  << exit(FatalError);
549  }
550 }
551 
552 
553 template<class Type, template<class> class PatchField, class GeoMesh>
556 types() const
557 {
558  const FieldField<PatchField, Type>& pff = *this;
559 
560  wordList Types(pff.size());
561 
562  forAll(pff, patchi)
563  {
564  Types[patchi] = pff[patchi].type();
565  }
566 
567  return Types;
568 }
569 
570 
571 template<class Type, template<class> class PatchField, class GeoMesh>
575 {
577  BoundaryInternalField(*this);
578 
579  forAll(BoundaryInternalField, patchi)
580  {
581  BoundaryInternalField[patchi] ==
582  this->operator[](patchi).patchInternalField();
583  }
584 
585  return BoundaryInternalField;
586 }
587 
588 
589 template<class Type, template<class> class PatchField, class GeoMesh>
592 interfaces() const
593 {
595 
596  forAll(interfaces, patchi)
597  {
598  if (isA<LduInterfaceField<Type> >(this->operator[](patchi)))
599  {
600  interfaces.set
601  (
602  patchi,
604  (
605  this->operator[](patchi)
606  )
607  );
608  }
609  }
610 
611  return interfaces;
612 }
613 
614 
615 template<class Type, template<class> class PatchField, class GeoMesh>
619 {
621 
622  forAll(interfaces, patchi)
623  {
624  if (isA<lduInterfaceField>(this->operator[](patchi)))
625  {
626  interfaces.set
627  (
628  patchi,
629  &refCast<const lduInterfaceField>
630  (
631  this->operator[](patchi)
632  )
633  );
634  }
635  }
636 
637  return interfaces;
638 }
639 
640 
641 template<class Type, template<class> class PatchField, class GeoMesh>
643 writeEntry(const word& keyword, Ostream& os) const
644 {
645  os << keyword << nl << token::BEGIN_BLOCK << incrIndent << nl;
646 
647  forAll(*this, patchi)
648  {
649  os << indent << this->operator[](patchi).patch().name() << nl
650  << indent << token::BEGIN_BLOCK << nl
651  << incrIndent << this->operator[](patchi) << decrIndent
652  << indent << token::END_BLOCK << endl;
653  }
654 
655  os << decrIndent << token::END_BLOCK << endl;
656 
657  // Check state of IOstream
658  os.check
659  (
660  "GeometricField<Type, PatchField, GeoMesh>::GeometricBoundaryField::"
661  "writeEntry(const word& keyword, Ostream& os) const"
662  );
663 }
664 
665 
666 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
667 
668 template<class Type, template<class> class PatchField, class GeoMesh>
670 operator=
671 (
673  GeometricBoundaryField& bf
674 )
675 {
677 }
678 
679 
680 template<class Type, template<class> class PatchField, class GeoMesh>
682 operator=
683 (
684  const FieldField<PatchField, Type>& ptff
685 )
686 {
688 }
689 
690 
691 template<class Type, template<class> class PatchField, class GeoMesh>
693 operator=
694 (
695  const Type& t
696 )
697 {
699 }
700 
701 
702 // Forced assignments
703 template<class Type, template<class> class PatchField, class GeoMesh>
705 operator==
706 (
708  GeometricBoundaryField& bf
709 )
710 {
711  forAll((*this), patchI)
712  {
713  this->operator[](patchI) == bf[patchI];
714  }
715 }
716 
717 
718 template<class Type, template<class> class PatchField, class GeoMesh>
720 operator==
721 (
722  const FieldField<PatchField, Type>& ptff
723 )
724 {
725  forAll((*this), patchI)
726  {
727  this->operator[](patchI) == ptff[patchI];
728  }
729 }
730 
731 
732 template<class Type, template<class> class PatchField, class GeoMesh>
734 operator==
735 (
736  const Type& t
737 )
738 {
739  forAll((*this), patchI)
740  {
741  this->operator[](patchI) == t;
742  }
743 }
744 
745 
746 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
747 
748 template<class Type, template<class> class PatchField, class GeoMesh>
749 Foam::Ostream& Foam::operator<<
750 (
751  Ostream& os,
753  GeometricBoundaryField& bf
754 )
755 {
756  os << static_cast<const FieldField<PatchField, Type>&>(bf);
757  return os;
758 }
759 
760 
761 // ************************************************************************* //
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:380
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
virtual bool isDict() const
Return true if this entry is a dictionary.
Definition: entry.H:153
void evaluate()
Evaluate boundary conditions.
label size() const
Return the number of elements in the PtrList.
Definition: PtrListI.H:32
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void writeEntry(const word &keyword, Ostream &os) const
Write boundary field as dictionary entry.
An abstract base class for implicitly-coupled interface fields e.g. processor and cyclic patch fields...
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual const dictionary & dict() const =0
Return dictionary if this entry is a dictionary.
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:221
friend Ostream & operator(Ostream &, const FieldField< PatchField, Type > &)
const_reverse_iterator rbegin() const
Definition: UILList.H:342
A class for handling words, derived from string.
Definition: word.H:59
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
bool isPattern() const
Should be treated as a match rather than a literal string.
Definition: keyTypeI.H:76
messageStream Info
LduInterfaceFieldPtrsList< Type > interfaces() const
Return a list of pointers for each patch field with only those.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:638
label size() const
Return number of elements in list.
Definition: DLListBaseI.H:77
Intrusive doubly-linked list.
Definition: IDLList.H:47
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
const keyType & keyword() const
Return keyword.
Definition: entry.H:120
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
static void waitRequests(const label start=0)
Wait until all requests (from start onwards) have finished.
Definition: UPstream.C:106
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
void readField(const DimensionedField< Type, GeoMesh > &field, const dictionary &dict)
Read the boundary field.
GeometricBoundaryField boundaryInternalField() const
Return BoundaryField of the cell values neighbouring.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:228
dictionary dict
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:235
IOerror FatalIOError
GeoMesh::BoundaryMesh BoundaryMesh
void updateCoeffs()
Update the boundary condition coefficients.
const dimensionedScalar e
Elementary charge.
Definition: doubleFloat.H:78
#define forAll(list, i)
Definition: UList.H:421
label patchi
errorManip< error > abort(error &err)
Definition: errorManip.H:131
tmp< DimensionedField< Type, GeoMesh > > clone() const
Clone.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
wordList types() const
Return a list of the patch types.
bool set(const label) const
Is element set.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
static label nRequests()
Get number of outstanding requests.
Definition: UPstream.C:96
bool found
A keyword and a list of tokens is an &#39;entry&#39;.
Definition: entry.H:65
points setSize(newPointi)
Generic GeometricField class.
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
error FatalError
static const NamedEnum< commsTypes, 3 > commsTypeNames
Definition: UPstream.H:71
lduInterfaceFieldPtrsList scalarInterfaces() const
Return a list of pointers for each patch field with only those.
void operator=(const FieldField< Field, Type > &)
Definition: FieldField.C:296
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
GeometricBoundaryField(const BoundaryMesh &)
Construct from a BoundaryMesh.
const const_reverse_iterator & rend() const
Definition: UILList.H:347
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:106
static commsTypes defaultCommsType
Default commsType.
Definition: UPstream.H:255
UEqn clear()
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
const PatchField< Type > & operator[](const label) const
Return element const reference.