SlicedGeometricField.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-2017 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 "SlicedGeometricField.H"
27 #include "processorFvPatch.H"
28 
29 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * * //
30 
31 template
32 <
33  class Type,
34  template<class> class PatchField,
35  template<class> class SlicedPatchField,
36  class GeoMesh
37 >
41 (
42  const Mesh& mesh,
43  const Field<Type>& completeField,
44  const bool preserveCouples,
45  const bool preserveProcessorOnly
46 )
47 {
48  tmp<FieldField<PatchField, Type>> tbf
49  (
50  new FieldField<PatchField, Type>(mesh.boundary().size())
51  );
52  FieldField<PatchField, Type>& bf = tbf.ref();
53 
54  forAll(mesh.boundary(), patchi)
55  {
56  if
57  (
58  preserveCouples
59  && mesh.boundary()[patchi].coupled()
60  && (
61  !preserveProcessorOnly
62  || isA<processorFvPatch>(mesh.boundary()[patchi])
63  )
64  )
65  {
66  // For coupled patched construct the correct patch field type
67  bf.set
68  (
69  patchi,
71  (
72  mesh.boundary()[patchi].type(),
73  mesh.boundary()[patchi],
74  *this
75  )
76  );
77 
78  // Initialize the values on the coupled patch to those of the slice
79  // of the given field.
80  // Note: these will usually be over-ridden by the boundary field
81  // evaluation e.g. in the case of processor and cyclic patches.
82  bf[patchi] = SlicedPatchField<Type>
83  (
84  mesh.boundary()[patchi],
85  DimensionedField<Type, GeoMesh>::null(),
86  completeField
87  );
88  }
89  else
90  {
91  bf.set
92  (
93  patchi,
94  new SlicedPatchField<Type>
95  (
96  mesh.boundary()[patchi],
97  DimensionedField<Type, GeoMesh>::null(),
98  completeField
99  )
100  );
101  }
102  }
103 
104  return tbf;
105 }
106 
107 
108 template
109 <
110  class Type,
111  template<class> class PatchField,
112  template<class> class SlicedPatchField,
113  class GeoMesh
114 >
118 (
119  const Mesh& mesh,
120  const FieldField<PatchField, Type>& bField,
121  const bool preserveCouples
122 )
123 {
124  tmp<FieldField<PatchField, Type>> tbf
125  (
126  new FieldField<PatchField, Type>(mesh.boundary().size())
127  );
128  FieldField<PatchField, Type>& bf = tbf.ref();
129 
130  forAll(mesh.boundary(), patchi)
131  {
132  if (preserveCouples && mesh.boundary()[patchi].coupled())
133  {
134  // For coupled patched construct the correct patch field type
135  bf.set
136  (
137  patchi,
139  (
140  mesh.boundary()[patchi].type(),
141  mesh.boundary()[patchi],
142  *this
143  )
144  );
145 
146  // Assign field
147  bf[patchi] == bField[patchi];
148  }
149  else
150  {
151  // Create unallocated copy of patch field
152  bf.set
153  (
154  patchi,
155  new SlicedPatchField<Type>
156  (
157  mesh.boundary()[patchi],
158  DimensionedField<Type, GeoMesh>::null()
159  )
160  );
161  bf[patchi].UList<Type>::shallowCopy(bField[patchi]);
162  }
163  }
164 
165  return tbf;
166 }
167 
168 
169 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
170 
171 template
172 <
173  class Type,
174  template<class> class PatchField,
175  template<class> class SlicedPatchField,
176  class GeoMesh
177 >
180 (
181  const IOobject& io,
182  const Mesh& mesh,
183  const dimensionSet& ds,
184  const Field<Type>& iField
185 )
186 :
188  (
189  io,
190  mesh,
191  ds,
192  Field<Type>()
193  )
194 {
195  // Set the internalField to the slice of the complete field
197  (
198  typename Field<Type>::subField(iField, GeoMesh::size(mesh))
199  );
200 }
201 
202 
203 template
204 <
205  class Type,
206  template<class> class PatchField,
207  template<class> class SlicedPatchField,
208  class GeoMesh
209 >
212 (
213  const IOobject& io,
214  const Mesh& mesh,
215  const dimensionSet& ds,
216  const Field<Type>& completeField,
217  const bool preserveCouples
218 )
219 :
221  (
222  io,
223  mesh,
224  ds,
225  Field<Type>(),
226  slicedBoundaryField(mesh, completeField, preserveCouples)
227  )
228 {
229  // Set the internalField to the slice of the complete field
231  (
232  typename Field<Type>::subField(completeField, GeoMesh::size(mesh))
233  );
234 
236 }
237 
238 
239 template
240 <
241  class Type,
242  template<class> class PatchField,
243  template<class> class SlicedPatchField,
244  class GeoMesh
245 >
248 (
249  const IOobject& io,
250  const Mesh& mesh,
251  const dimensionSet& ds,
252  const Field<Type>& completeIField,
253  const Field<Type>& completeBField,
254  const bool preserveCouples,
255  const bool preserveProcessorOnly
256 )
257 :
259  (
260  io,
261  mesh,
262  ds,
263  Field<Type>(),
264  slicedBoundaryField
265  (
266  mesh,
267  completeBField,
268  preserveCouples,
269  preserveProcessorOnly
270  )
271  )
272 {
273  // Set the internalField to the slice of the complete field
275  (
276  typename Field<Type>::subField(completeIField, GeoMesh::size(mesh))
277  );
278 
280 }
281 
282 
283 template
284 <
285  class Type,
286  template<class> class PatchField,
287  template<class> class SlicedPatchField,
288  class GeoMesh
289 >
292 (
293  const IOobject& io,
295  const bool preserveCouples
296 )
297 :
299  (
300  io,
301  gf.mesh(),
302  gf.dimensions(),
303  Field<Type>(),
304  slicedBoundaryField(gf.mesh(), gf.boundaryField(), preserveCouples)
305  )
306 {
307  // Set the internalField to the supplied internal field
309 
311 }
312 
313 
314 template
315 <
316  class Type,
317  template<class> class PatchField,
318  template<class> class SlicedPatchField,
319  class GeoMesh
320 >
323 (
325 )
326 :
328  (
329  gf,
330  gf.mesh(),
331  gf.dimensions(),
332  Field<Type>(),
333  slicedBoundaryField(gf.mesh(), gf.boundaryField(), true)
334  )
335 {
336  // Set the internalField to the supplied internal field
338 }
339 
340 
341 template
342 <
343  class Type,
344  template<class> class PatchField,
345  template<class> class SlicedPatchField,
346  class GeoMesh
347 >
348 Foam::tmp
349 <
351 >
353 clone() const
354 {
355  return tmp
356  <
358  >
359  (
361  (
362  *this
363  )
364  );
365 }
366 
367 
368 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
369 
370 template
371 <
372  class Type,
373  template<class> class PatchField,
374  template<class> class SlicedPatchField,
375  class GeoMesh
376 >
379 {
380  // Set the internalField storage pointer to nullptr before its destruction
381  // to protect the field it a slice of.
383 }
384 
385 
386 template
387 <
388  class Type,
389  template<class> class PatchField,
390  template<class> class SlicedPatchField,
391  class GeoMesh
392 >
395 {
396  // Set the internalField storage pointer to nullptr before its destruction
397  // to protect the field it a slice of.
399 }
400 
401 
402 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
403 
404 template
405 <
406  class Type,
407  template<class> class PatchField,
408  template<class> class SlicedPatchField,
409  class GeoMesh
410 >
413 {
415 }
416 
417 
418 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
Pre-declare related SubField type.
Definition: Field.H:61
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
const dimensionSet & dimensions() const
Return dimensions.
Dimension set for the base types.
Definition: dimensionSet.H:120
dynamicFvMesh & mesh
Pre-declare SubField and related Field type.
Definition: Field.H:57
GeoMesh::Mesh Mesh
Type of mesh on which this DimensionedField is instantiated.
Internal(const IOobject &, const Mesh &, const dimensionSet &, const Field< Type > &iField)
Construct from components and field to slice.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:61
const Mesh & mesh() const
Return mesh.
tmp< SlicedGeometricField< Type, PatchField, SlicedPatchField, GeoMesh > > clone() const
Clone.
void correctBoundaryConditions()
Correct boundary field.
label patchi
U correctBoundaryConditions()
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
SlicedGeometricField(const IOobject &, const Mesh &, const dimensionSet &, const Field< Type > &completeField, const bool preserveCouples=true)
Construct from components and field to slice.
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
A class for managing temporary objects.
Definition: PtrList.H:53
Specialization of GeometricField which holds slices of given complete fields in a form that they act ...
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92