genericFvPatchField.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-2022 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 "genericFvPatchField.H"
27 #include "fvPatchFieldMapper.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Type>
33 (
34  const fvPatch& p,
36 )
37 :
39 {
41  << "Trying to construct an genericFvPatchField on patch "
42  << this->patch().name()
43  << " of field " << this->internalField().name()
44  << abort(FatalError);
45 }
46 
47 
48 template<class Type>
50 (
51  const fvPatch& p,
53  const dictionary& dict
54 )
55 :
56  genericPatchField(dict.lookup("type")),
57  calculatedFvPatchField<Type>(p, iF, dict),
58  dict_(dict)
59 {
60  if (!dict.found("value"))
61  {
63  (
64  dict
65  ) << "\n Cannot find 'value' entry"
66  << " on patch " << this->patch().name()
67  << " of field " << this->internalField().name()
68  << " in file " << this->internalField().objectPath()
69  << nl
70  << " which is required to set the"
71  " values of the generic patch field." << nl
72  << " (Actual type " << actualTypeName() << ")" << nl
73  << "\n Please add the 'value' entry to the write function "
74  "of the user-defined boundary-condition\n"
75  << exit(FatalIOError);
76  }
77 
78  forAllConstIter(dictionary, dict_, iter)
79  {
80  if (iter().keyword() != "type" && iter().keyword() != "value")
81  {
82  if
83  (
84  iter().isStream()
85  && iter().stream().size()
86  )
87  {
88  ITstream& is = iter().stream();
89 
90  // Read first token
91  token firstToken(is);
92 
93  if
94  (
95  firstToken.isWord()
96  && firstToken.wordToken() == "nonuniform"
97  )
98  {
99  token fieldToken(is);
100 
101  if (!fieldToken.isCompound())
102  {
103  if
104  (
105  fieldToken.isLabel()
106  && fieldToken.labelToken() == 0
107  )
108  {
109  scalarFields_.insert
110  (
111  iter().keyword(),
112  new scalarField(0)
113  );
114  }
115  else
116  {
118  (
119  dict
120  ) << "\n token following 'nonuniform' "
121  "is not a compound"
122  << "\n on patch " << this->patch().name()
123  << " of field "
124  << this->internalField().name()
125  << " in file "
126  << this->internalField().objectPath()
127  << exit(FatalIOError);
128  }
129  }
130  else if
131  (
132  fieldToken.compoundToken().type()
133  == token::Compound<List<scalar>>::typeName
134  )
135  {
136  scalarField* fPtr = new scalarField;
137  fPtr->transfer
138  (
140  (
141  fieldToken.transferCompoundToken(is)
142  )
143  );
144 
145  if (fPtr->size() != this->size())
146  {
148  (
149  dict
150  ) << "\n size of field " << iter().keyword()
151  << " (" << fPtr->size() << ')'
152  << " is not the same size as the patch ("
153  << this->size() << ')'
154  << "\n on patch " << this->patch().name()
155  << " of field "
156  << this->internalField().name()
157  << " in file "
158  << this->internalField().objectPath()
159  << exit(FatalIOError);
160  }
161 
162  scalarFields_.insert(iter().keyword(), fPtr);
163  }
164  else if
165  (
166  fieldToken.compoundToken().type()
167  == token::Compound<List<vector>>::typeName
168  )
169  {
170  vectorField* fPtr = new vectorField;
171  fPtr->transfer
172  (
174  (
175  fieldToken.transferCompoundToken(is)
176  )
177  );
178 
179  if (fPtr->size() != this->size())
180  {
182  (
183  dict
184  ) << "\n size of field " << iter().keyword()
185  << " (" << fPtr->size() << ')'
186  << " is not the same size as the patch ("
187  << this->size() << ')'
188  << "\n on patch " << this->patch().name()
189  << " of field "
190  << this->internalField().name()
191  << " in file "
192  << this->internalField().objectPath()
193  << exit(FatalIOError);
194  }
195 
196  vectorFields_.insert(iter().keyword(), fPtr);
197  }
198  else if
199  (
200  fieldToken.compoundToken().type()
202  )
203  {
205  fPtr->transfer
206  (
208  <
210  >
211  (
212  fieldToken.transferCompoundToken(is)
213  )
214  );
215 
216  if (fPtr->size() != this->size())
217  {
219  (
220  dict
221  ) << "\n size of field " << iter().keyword()
222  << " (" << fPtr->size() << ')'
223  << " is not the same size as the patch ("
224  << this->size() << ')'
225  << "\n on patch " << this->patch().name()
226  << " of field "
227  << this->internalField().name()
228  << " in file "
229  << this->internalField().objectPath()
230  << exit(FatalIOError);
231  }
232 
233  sphericalTensorFields_.insert(iter().keyword(), fPtr);
234  }
235  else if
236  (
237  fieldToken.compoundToken().type()
239  )
240  {
241  symmTensorField* fPtr = new symmTensorField;
242  fPtr->transfer
243  (
245  <
247  >
248  (
249  fieldToken.transferCompoundToken(is)
250  )
251  );
252 
253  if (fPtr->size() != this->size())
254  {
256  (
257  dict
258  ) << "\n size of field " << iter().keyword()
259  << " (" << fPtr->size() << ')'
260  << " is not the same size as the patch ("
261  << this->size() << ')'
262  << "\n on patch " << this->patch().name()
263  << " of field "
264  << this->internalField().name()
265  << " in file "
266  << this->internalField().objectPath()
267  << exit(FatalIOError);
268  }
269 
270  symmTensorFields_.insert(iter().keyword(), fPtr);
271  }
272  else if
273  (
274  fieldToken.compoundToken().type()
275  == token::Compound<List<tensor>>::typeName
276  )
277  {
278  tensorField* fPtr = new tensorField;
279  fPtr->transfer
280  (
282  (
283  fieldToken.transferCompoundToken(is)
284  )
285  );
286 
287  if (fPtr->size() != this->size())
288  {
290  (
291  dict
292  ) << "\n size of field " << iter().keyword()
293  << " (" << fPtr->size() << ')'
294  << " is not the same size as the patch ("
295  << this->size() << ')'
296  << "\n on patch " << this->patch().name()
297  << " of field "
298  << this->internalField().name()
299  << " in file "
300  << this->internalField().objectPath()
301  << exit(FatalIOError);
302  }
303 
304  tensorFields_.insert(iter().keyword(), fPtr);
305  }
306  else
307  {
309  (
310  dict
311  ) << "\n compound " << fieldToken.compoundToken()
312  << " not supported"
313  << "\n on patch " << this->patch().name()
314  << " of field "
315  << this->internalField().name()
316  << " in file "
317  << this->internalField().objectPath()
318  << exit(FatalIOError);
319  }
320  }
321  else if
322  (
323  firstToken.isWord()
324  && firstToken.wordToken() == "uniform"
325  )
326  {
327  token fieldToken(is);
328 
329  if (!fieldToken.isPunctuation())
330  {
331  scalarFields_.insert
332  (
333  iter().keyword(),
334  new scalarField
335  (
336  this->size(),
337  fieldToken.number()
338  )
339  );
340  }
341  else
342  {
343  // Read as scalarList.
344  is.putBack(fieldToken);
345 
346  scalarList l(is);
347 
348  if (l.size() == vector::nComponents)
349  {
350  vector vs(l[0], l[1], l[2]);
351 
352  vectorFields_.insert
353  (
354  iter().keyword(),
355  new vectorField(this->size(), vs)
356  );
357  }
358  else if (l.size() == sphericalTensor::nComponents)
359  {
360  sphericalTensor vs(l[0]);
361 
362  sphericalTensorFields_.insert
363  (
364  iter().keyword(),
365  new sphericalTensorField(this->size(), vs)
366  );
367  }
368  else if (l.size() == symmTensor::nComponents)
369  {
370  symmTensor vs(l[0], l[1], l[2], l[3], l[4], l[5]);
371 
372  symmTensorFields_.insert
373  (
374  iter().keyword(),
375  new symmTensorField(this->size(), vs)
376  );
377  }
378  else if (l.size() == tensor::nComponents)
379  {
380  tensor vs
381  (
382  l[0], l[1], l[2],
383  l[3], l[4], l[5],
384  l[6], l[7], l[8]
385  );
386 
387  tensorFields_.insert
388  (
389  iter().keyword(),
390  new tensorField(this->size(), vs)
391  );
392  }
393  else
394  {
396  (
397  dict
398  ) << "\n unrecognised native type " << l
399  << "\n on patch " << this->patch().name()
400  << " of field "
401  << this->internalField().name()
402  << " in file "
403  << this->internalField().objectPath()
404  << exit(FatalIOError);
405  }
406  }
407  }
408  }
409  }
410  }
411 }
412 
413 
414 template<class Type>
416 (
417  const genericFvPatchField<Type>& ptf,
418  const fvPatch& p,
420  const fvPatchFieldMapper& mapper
421 )
422 :
423  genericPatchField(ptf),
424  calculatedFvPatchField<Type>(ptf, p, iF, mapper),
425  dict_(ptf.dict_)
426 {
428  (
430  ptf.scalarFields_,
431  iter
432  )
433  {
434  scalarFields_.insert
435  (
436  iter.key(),
437  mapper(*iter()).ptr()
438  );
439  }
440 
442  (
444  ptf.vectorFields_,
445  iter
446  )
447  {
448  vectorFields_.insert
449  (
450  iter.key(),
451  mapper(*iter()).ptr()
452  );
453  }
454 
456  (
458  ptf.sphericalTensorFields_,
459  iter
460  )
461  {
462  sphericalTensorFields_.insert
463  (
464  iter.key(),
465  mapper(*iter()).ptr()
466  );
467  }
468 
470  (
472  ptf.symmTensorFields_,
473  iter
474  )
475  {
476  symmTensorFields_.insert
477  (
478  iter.key(),
479  mapper(*iter()).ptr()
480  );
481  }
482 
484  (
486  ptf.tensorFields_,
487  iter
488  )
489  {
490  tensorFields_.insert
491  (
492  iter.key(),
493  mapper(*iter()).ptr()
494  );
495  }
496 }
497 
498 
499 template<class Type>
501 (
502  const genericFvPatchField<Type>& ptf,
504 )
505 :
506  genericPatchField(ptf),
508  dict_(ptf.dict_),
509  scalarFields_(ptf.scalarFields_),
510  vectorFields_(ptf.vectorFields_),
511  sphericalTensorFields_(ptf.sphericalTensorFields_),
512  symmTensorFields_(ptf.symmTensorFields_),
513  tensorFields_(ptf.tensorFields_)
514 {}
515 
516 
517 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
518 
519 template<class Type>
521 (
522  const fvPatchFieldMapper& m
523 )
524 {
526 
527  forAllIter
528  (
530  scalarFields_,
531  iter
532  )
533  {
534  m(*iter(), *iter());
535  }
536 
537  forAllIter
538  (
540  vectorFields_,
541  iter
542  )
543  {
544  m(*iter(), *iter());
545  }
546 
547  forAllIter
548  (
550  sphericalTensorFields_,
551  iter
552  )
553  {
554  m(*iter(), *iter());
555  }
556 
557  forAllIter
558  (
560  symmTensorFields_,
561  iter
562  )
563  {
564  m(*iter(), *iter());
565  }
566 
567  forAllIter
568  (
570  tensorFields_,
571  iter
572  )
573  {
574  m(*iter(), *iter());
575  }
576 }
577 
578 
579 template<class Type>
581 (
582  const fvPatchField<Type>& ptf,
583  const labelList& addr
584 )
585 {
587 
588  const genericFvPatchField<Type>& dptf =
589  refCast<const genericFvPatchField<Type>>(ptf);
590 
591  forAllIter
592  (
594  scalarFields_,
595  iter
596  )
597  {
599  dptf.scalarFields_.find(iter.key());
600 
601  if (dptfIter != dptf.scalarFields_.end())
602  {
603  iter()->rmap(*dptfIter(), addr);
604  }
605  }
606 
607  forAllIter
608  (
610  vectorFields_,
611  iter
612  )
613  {
615  dptf.vectorFields_.find(iter.key());
616 
617  if (dptfIter != dptf.vectorFields_.end())
618  {
619  iter()->rmap(*dptfIter(), addr);
620  }
621  }
622 
623  forAllIter
624  (
626  sphericalTensorFields_,
627  iter
628  )
629  {
631  dptf.sphericalTensorFields_.find(iter.key());
632 
633  if (dptfIter != dptf.sphericalTensorFields_.end())
634  {
635  iter()->rmap(*dptfIter(), addr);
636  }
637  }
638 
639  forAllIter
640  (
642  symmTensorFields_,
643  iter
644  )
645  {
647  dptf.symmTensorFields_.find(iter.key());
648 
649  if (dptfIter != dptf.symmTensorFields_.end())
650  {
651  iter()->rmap(*dptfIter(), addr);
652  }
653  }
654 
655  forAllIter
656  (
658  tensorFields_,
659  iter
660  )
661  {
663  dptf.tensorFields_.find(iter.key());
664 
665  if (dptfIter != dptf.tensorFields_.end())
666  {
667  iter()->rmap(*dptfIter(), addr);
668  }
669  }
670 }
671 
672 
673 template<class Type>
675 (
676  const fvPatchField<Type>& ptf
677 )
678 {
680 
681  const genericFvPatchField<Type>& dptf =
682  refCast<const genericFvPatchField<Type>>(ptf);
683 
684  forAllIter
685  (
687  scalarFields_,
688  iter
689  )
690  {
692  dptf.scalarFields_.find(iter.key());
693 
694  if (dptfIter != dptf.scalarFields_.end())
695  {
696  iter()->reset(*dptfIter());
697  }
698  }
699 
700  forAllIter
701  (
703  vectorFields_,
704  iter
705  )
706  {
708  dptf.vectorFields_.find(iter.key());
709 
710  if (dptfIter != dptf.vectorFields_.end())
711  {
712  iter()->reset(*dptfIter());
713  }
714  }
715 
716  forAllIter
717  (
719  sphericalTensorFields_,
720  iter
721  )
722  {
724  dptf.sphericalTensorFields_.find(iter.key());
725 
726  if (dptfIter != dptf.sphericalTensorFields_.end())
727  {
728  iter()->reset(*dptfIter());
729  }
730  }
731 
732  forAllIter
733  (
735  symmTensorFields_,
736  iter
737  )
738  {
740  dptf.symmTensorFields_.find(iter.key());
741 
742  if (dptfIter != dptf.symmTensorFields_.end())
743  {
744  iter()->reset(*dptfIter());
745  }
746  }
747 
748  forAllIter
749  (
751  tensorFields_,
752  iter
753  )
754  {
756  dptf.tensorFields_.find(iter.key());
757 
758  if (dptfIter != dptf.tensorFields_.end())
759  {
760  iter()->reset(*dptfIter());
761  }
762  }
763 }
764 
765 
766 template<class Type>
769 (
770  const tmp<scalarField>&
771 ) const
772 {
774  << "cannot be called for a genericFvPatchField"
775  " (actual type " << actualTypeName() << ")"
776  << "\n on patch " << this->patch().name()
777  << " of field " << this->internalField().name()
778  << " in file " << this->internalField().objectPath()
779  << "\n You are probably trying to solve for a field with a "
780  "generic boundary condition."
781  << abort(FatalError);
782 
783  return *this;
784 }
785 
786 
787 template<class Type>
790 (
791  const tmp<scalarField>&
792 ) const
793 {
795  << "cannot be called for a genericFvPatchField"
796  " (actual type " << actualTypeName() << ")"
797  << "\n on patch " << this->patch().name()
798  << " of field " << this->internalField().name()
799  << " in file " << this->internalField().objectPath()
800  << "\n You are probably trying to solve for a field with a "
801  "generic boundary condition."
802  << abort(FatalError);
803 
804  return *this;
805 }
806 
807 
808 template<class Type>
811 {
813  << "cannot be called for a genericFvPatchField"
814  " (actual type " << actualTypeName() << ")"
815  << "\n on patch " << this->patch().name()
816  << " of field " << this->internalField().name()
817  << " in file " << this->internalField().objectPath()
818  << "\n You are probably trying to solve for a field with a "
819  "generic boundary condition."
820  << abort(FatalError);
821 
822  return *this;
823 }
824 
825 template<class Type>
828 {
830  << "cannot be called for a genericFvPatchField"
831  " (actual type " << actualTypeName() << ")"
832  << "\n on patch " << this->patch().name()
833  << " of field " << this->internalField().name()
834  << " in file " << this->internalField().objectPath()
835  << "\n You are probably trying to solve for a field with a "
836  "generic boundary condition."
837  << abort(FatalError);
838 
839  return *this;
840 }
841 
842 
843 template<class Type>
845 {
846  writeEntry(os, "type", actualTypeName());
847 
848  forAllConstIter(dictionary, dict_, iter)
849  {
850  if (iter().keyword() != "type" && iter().keyword() != "value")
851  {
852  if
853  (
854  iter().isStream()
855  && iter().stream().size()
856  && iter().stream()[0].isWord()
857  && iter().stream()[0].wordToken() == "nonuniform"
858  )
859  {
860  if (scalarFields_.found(iter().keyword()))
861  {
862  writeEntry
863  (
864  os,
865  iter().keyword(),
866  *scalarFields_.find(iter().keyword())()
867  );
868  }
869  else if (vectorFields_.found(iter().keyword()))
870  {
871  writeEntry
872  (
873  os,
874  iter().keyword(),
875  *vectorFields_.find(iter().keyword())()
876  );
877  }
878  else if (sphericalTensorFields_.found(iter().keyword()))
879  {
880  writeEntry
881  (
882  os,
883  iter().keyword(),
884  *sphericalTensorFields_.find(iter().keyword())()
885  );
886  }
887  else if (symmTensorFields_.found(iter().keyword()))
888  {
889  writeEntry
890  (
891  os,
892  iter().keyword(),
893  *symmTensorFields_.find(iter().keyword())()
894  );
895  }
896  else if (tensorFields_.found(iter().keyword()))
897  {
898  writeEntry
899  (
900  os,
901  iter().keyword(),
902  *tensorFields_.find(iter().keyword())()
903  );
904  }
905  }
906  else
907  {
908  iter().write(os);
909  }
910  }
911  }
912 
913  writeEntry(os, "value", *this);
914 }
915 
916 
917 // ************************************************************************* //
bool isLabel() const
Definition: tokenI.H:392
tmp< Field< Type > > gradientInternalCoeffs() const
Return the matrix diagonal coefficients corresponding to the.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:663
bool isWord() const
Definition: tokenI.H:261
Field< sphericalTensor > sphericalTensorField
Specialisation of Field<T> for sphericalTensor.
compound & transferCompoundToken(const Istream &is)
Definition: token.C:93
Field< symmTensor > symmTensorField
Specialisation of Field<T> for symmTensor.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
To & dynamicCast(From &r)
Reference type cast template function,.
Definition: typeInfo.H:85
An STL-conforming const_iterator.
Definition: HashTable.H:481
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:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
This boundary condition provides a generic version of the calculated condition, useful as a fallback ...
scalar number() const
Definition: tokenI.H:503
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:459
const word & wordToken() const
Definition: tokenI.H:266
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
A token holds items read from Istream.
Definition: token.H:72
void putBack(const token &)
Put back token.
Definition: Istream.C:30
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:63
virtual void write(Ostream &) const
Write.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
volVectorField vectorField(fieldObject, mesh)
A HashTable specialisation for hashing pointers.
Definition: HashPtrTable.H:50
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
bool isCompound() const
Definition: tokenI.H:520
Templated 3D SphericalTensor derived from VectorSpace adding construction from 1 component, element access using th ii() member function and the inner-product (dot-product) and outer-product operators.
A templated class for holding compound tokens.
Definition: token.H:209
virtual tmp< Field< Type > > valueBoundaryCoeffs(const tmp< scalarField > &) const
Return the matrix source coefficients corresponding to the.
virtual void reset(const fvPatchField< Type > &)
Reset the fvPatchField to the given fvPatchField.
Foam::fvPatchFieldMapper.
tmp< Field< Type > > gradientBoundaryCoeffs() const
Return the matrix source coefficients corresponding to the.
genericFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
This boundary condition is not designed to be evaluated; it is assumed that the value is assigned via...
const compound & compoundToken() const
Definition: tokenI.H:525
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:54
volScalarField scalarField(fieldObject, mesh)
static const char nl
Definition: Ostream.H:260
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:318
label labelToken() const
Definition: tokenI.H:397
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual tmp< Field< Type > > valueInternalCoeffs(const tmp< scalarField > &) const
Return the matrix diagonal coefficients corresponding to the.
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342
bool isPunctuation() const
Definition: tokenI.H:243
Input token stream.
Definition: ITstream.H:49
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864
IOerror FatalIOError