FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 
27 #include <OpenFOAM/commSchedule.H>
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
32 template<class Type, template<class> class PatchField, class GeoMesh>
35 (
36  const BoundaryMesh& bmesh,
37  const DimensionedField<Type, GeoMesh>& field,
38  const word& patchFieldType
39 )
40 :
41  FieldField<PatchField, Type>(bmesh.size()),
42  bmesh_(bmesh)
43 {
44  if (debug)
45  {
46  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
47  "GeometricBoundaryField::"
48  "GeometricBoundaryField(const BoundaryMesh&, "
49  "const Field<Type>&, const word&)"
50  << endl;
51  }
52 
53  forAll(bmesh_, patchi)
54  {
55  set
56  (
57  patchi,
58  PatchField<Type>::New
59  (
60  patchFieldType,
61  bmesh_[patchi],
62  field
63  )
64  );
65  }
66 }
67 
68 
69 template<class Type, template<class> class PatchField, class GeoMesh>
72 (
73  const BoundaryMesh& bmesh,
74  const DimensionedField<Type, GeoMesh>& field,
75  const wordList& patchFieldTypes
76 )
77 :
78  FieldField<PatchField, Type>(bmesh.size()),
79  bmesh_(bmesh)
80 {
81  if (debug)
82  {
83  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
84  "GeometricBoundaryField::"
85  "GeometricBoundaryField(const BoundaryMesh&, "
86  "const Field<Type>&, const wordList&)"
87  << endl;
88  }
89 
90  if (patchFieldTypes.size() != this->size())
91  {
93  (
94  "GeometricField<Type, PatchField, GeoMesh>::"
95  "GeometricBoundaryField::"
96  "GeometricBoundaryField(const BoundaryMesh&, "
97  "const Field<Type>&, const wordList&)"
98  ) << "Incorrect number of patch type specifications given" << nl
99  << " Number of patches in mesh = " << bmesh.size()
100  << " number of patch type specifications = "
101  << patchFieldTypes.size()
102  << abort(FatalError);
103  }
104 
105  forAll(bmesh_, patchi)
106  {
107  set
108  (
109  patchi,
110  PatchField<Type>::New
111  (
112  patchFieldTypes[patchi],
113  bmesh_[patchi],
114  field
115  )
116  );
117  }
118 }
119 
120 
121 template<class Type, template<class> class PatchField, class GeoMesh>
124 (
125  const BoundaryMesh& bmesh,
126  const DimensionedField<Type, GeoMesh>& field,
127  const PtrList<PatchField<Type> >& ptfl
128 )
129 :
130  FieldField<PatchField, Type>(bmesh.size()),
131  bmesh_(bmesh)
132 {
133  if (debug)
134  {
135  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
136  "GeometricBoundaryField::"
137  "GeometricBoundaryField(const BoundaryMesh&, "
138  "const Field<Type>&, const PatchField<Type>List&)"
139  << endl;
140  }
141 
142  forAll(bmesh_, patchi)
143  {
144  set(patchi, ptfl[patchi].clone(field));
145  }
146 }
147 
148 
149 template<class Type, template<class> class PatchField, class GeoMesh>
152 (
153  const DimensionedField<Type, GeoMesh>& field,
154  const typename GeometricField<Type, PatchField, GeoMesh>::
155  GeometricBoundaryField& btf
156 )
157 :
158  FieldField<PatchField, Type>(btf.size()),
159  bmesh_(btf.bmesh_)
160 {
161  if (debug)
162  {
163  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
164  "GeometricBoundaryField::"
165  "GeometricBoundaryField(const GeometricBoundaryField<Type, "
166  "PatchField, BoundaryMesh>&)"
167  << endl;
168  }
169 
170  forAll(bmesh_, patchi)
171  {
172  set(patchi, btf[patchi].clone(field));
173  }
174 }
175 
176 
177 // Construct as copy
178 // Dangerous because Field may be set to a field which gets deleted.
179 // Need new type of GeometricBoundaryField, one which IS part of a geometric
180 // field for which snGrad etc. may be called and a free standing
181 // GeometricBoundaryField for which such operations are unavailable.
182 template<class Type, template<class> class PatchField, class GeoMesh>
185 (
186  const typename GeometricField<Type, PatchField, GeoMesh>::
187  GeometricBoundaryField& btf
188 )
189 :
190  FieldField<PatchField, Type>(btf),
191  bmesh_(btf.bmesh_)
192 {
193  if (debug)
194  {
195  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
196  "GeometricBoundaryField::"
197  "GeometricBoundaryField(const GeometricBoundaryField<Type, "
198  "PatchField, BoundaryMesh>&)"
199  << endl;
200  }
201 }
202 
203 
204 template<class Type, template<class> class PatchField, class GeoMesh>
207 (
208  const BoundaryMesh& bmesh,
209  const DimensionedField<Type, GeoMesh>& field,
210  const dictionary& dict
211 )
212 :
213  FieldField<PatchField, Type>(bmesh.size()),
214  bmesh_(bmesh)
215 {
216  if (debug)
217  {
218  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
219  "GeometricBoundaryField::"
220  "GeometricBoundaryField"
221  "(const BoundaryMesh&, const Field<Type>&, const dictionary&)"
222  << endl;
223  }
224 
225  forAll(bmesh_, patchi)
226  {
227  if (bmesh_[patchi].type() != emptyPolyPatch::typeName)
228  {
229  set
230  (
231  patchi,
232  PatchField<Type>::New
233  (
234  bmesh_[patchi],
235  field,
236  dict.subDict(bmesh_[patchi].name())
237  )
238  );
239  }
240  else
241  {
242  set
243  (
244  patchi,
245  PatchField<Type>::New
246  (
247  emptyPolyPatch::typeName,
248  bmesh_[patchi],
249  field
250  )
251  );
252  }
253  }
254 }
255 
256 
257 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
258 
259 template<class Type, template<class> class PatchField, class GeoMesh>
262 {
263  if (debug)
264  {
265  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
266  "GeometricBoundaryField::"
267  "updateCoeffs()" << endl;
268  }
269 
270  forAll(*this, patchi)
271  {
272  this->operator[](patchi).updateCoeffs();
273  }
274 }
275 
276 
277 template<class Type, template<class> class PatchField, class GeoMesh>
280 {
281  if (debug)
282  {
283  Info<< "GeometricField<Type, PatchField, GeoMesh>::"
284  "GeometricBoundaryField::"
285  "evaluate()" << endl;
286  }
287 
288  if
289  (
292  )
293  {
294  forAll(*this, patchi)
295  {
296  this->operator[](patchi).initEvaluate(Pstream::defaultCommsType);
297  }
298 
299  // Block for any outstanding requests
301  {
304  }
305 
306  forAll(*this, patchi)
307  {
308  this->operator[](patchi).evaluate(Pstream::defaultCommsType);
309  }
310  }
312  {
313  const lduSchedule& patchSchedule =
314  bmesh_.mesh().globalData().patchSchedule();
315 
316  forAll(patchSchedule, patchEvali)
317  {
318  if (patchSchedule[patchEvali].init)
319  {
320  this->operator[](patchSchedule[patchEvali].patch)
321  .initEvaluate(Pstream::scheduled);
322  }
323  else
324  {
325  this->operator[](patchSchedule[patchEvali].patch)
326  .evaluate(Pstream::scheduled);
327  }
328  }
329  }
330  else
331  {
332  FatalErrorIn("GeometricBoundaryField::evaluate()")
333  << "Unsuported communications type "
335  << exit(FatalError);
336  }
337 }
338 
339 
340 template<class Type, template<class> class PatchField, class GeoMesh>
343 types() const
344 {
345  const FieldField<PatchField, Type>& pff = *this;
346 
347  wordList Types(pff.size());
348 
349  forAll(pff, patchi)
350  {
351  Types[patchi] = pff[patchi].type();
352  }
353 
354  return Types;
355 }
356 
357 
358 template<class Type, template<class> class PatchField, class GeoMesh>
362 {
364  BoundaryInternalField(*this);
365 
366  forAll(BoundaryInternalField, patchi)
367  {
368  BoundaryInternalField[patchi] ==
369  this->operator[](patchi).patchInternalField();
370  }
371 
372  return BoundaryInternalField;
373 }
374 
375 
376 template<class Type, template<class> class PatchField, class GeoMesh>
379 interfaces() const
380 {
381  lduInterfaceFieldPtrsList interfaces(this->size());
382 
383  forAll (interfaces, patchi)
384  {
385  if (isA<lduInterfaceField>(this->operator[](patchi)))
386  {
387  interfaces.set
388  (
389  patchi,
390  &refCast<const lduInterfaceField>(this->operator[](patchi))
391  );
392  }
393  }
394 
395  return interfaces;
396 }
397 
398 
399 template<class Type, template<class> class PatchField, class GeoMesh>
401 writeEntry(const word& keyword, Ostream& os) const
402 {
403  os << keyword << nl << token::BEGIN_BLOCK << incrIndent << nl;
404 
405  forAll(*this, patchi)
406  {
407  os << indent << this->operator[](patchi).patch().name() << nl
408  << indent << token::BEGIN_BLOCK << nl
409  << incrIndent << this->operator[](patchi) << decrIndent
410  << indent << token::END_BLOCK << endl;
411  }
412 
413  os << decrIndent << token::END_BLOCK << endl;
414 
415  // Check state of IOstream
416  os.check
417  (
418  "GeometricField<Type, PatchField, GeoMesh>::GeometricBoundaryField::"
419  "writeEntry(const word& keyword, Ostream& os) const"
420  );
421 }
422 
423 
424 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
425 
426 template<class Type, template<class> class PatchField, class GeoMesh>
428 operator=
429 (
432 )
433 {
435 }
436 
437 
438 template<class Type, template<class> class PatchField, class GeoMesh>
440 operator=
441 (
442  const FieldField<PatchField, Type>& ptff
443 )
444 {
446 }
447 
448 
449 template<class Type, template<class> class PatchField, class GeoMesh>
451 operator=
452 (
453  const Type& t
454 )
455 {
457 }
458 
459 
460 // Forced assignments
461 template<class Type, template<class> class PatchField, class GeoMesh>
463 operator==
464 (
467 )
468 {
469  forAll((*this), patchI)
470  {
471  this->operator[](patchI) == bf[patchI];
472  }
473 }
474 
475 
476 template<class Type, template<class> class PatchField, class GeoMesh>
478 operator==
479 (
480  const FieldField<PatchField, Type>& ptff
481 )
482 {
483  forAll((*this), patchI)
484  {
485  this->operator[](patchI) == ptff[patchI];
486  }
487 }
488 
489 
490 template<class Type, template<class> class PatchField, class GeoMesh>
492 operator==
493 (
494  const Type& t
495 )
496 {
497  forAll((*this), patchI)
498  {
499  this->operator[](patchI) == t;
500  }
501 }
502 
503 
504 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
505 
506 template<class Type, template<class> class PatchField, class GeoMesh>
507 Foam::Ostream& Foam::operator<<
508 (
509  Ostream& os,
511  GeometricBoundaryField& bf
512 )
513 {
514  os << static_cast<const FieldField<PatchField, Type>&>(bf);
515  return os;
516 }
517 
518 
519 // ************************ vim: set sw=4 sts=4 et: ************************ //