fvFieldDecomposerDecomposeFields.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 "fvFieldDecomposer.H"
27 #include "processorFvPatchField.H"
28 #include "processorFvsPatchField.H"
31 #include "emptyFvPatchFields.H"
32 
33 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
34 
35 template<class Type>
38 (
40  const bool allowUnknownPatchFields
41 ) const
42 {
43  // 1. Create the complete field with dummy patch fields
44  PtrList<fvPatchField<Type> > patchFields(boundaryAddressing_.size());
45 
46  forAll(boundaryAddressing_, patchi)
47  {
48  patchFields.set
49  (
50  patchi,
52  (
54  procMesh_.boundary()[patchi],
56  )
57  );
58  }
59 
60  // Create the field for the processor
62  (
64  (
65  IOobject
66  (
67  field.name(),
68  procMesh_.time().timeName(),
69  procMesh_,
70  IOobject::NO_READ,
71  IOobject::NO_WRITE
72  ),
73  procMesh_,
74  field.dimensions(),
75  Field<Type>(field.internalField(), cellAddressing_),
76  patchFields
77  )
78  );
80 
81 
82  // 2. Change the fvPatchFields to the correct type using a mapper
83  // constructor (with reference to the now correct internal field)
84 
86  GeometricBoundaryField& bf = resF.boundaryField();
87 
88  forAll(bf, patchi)
89  {
90  if (patchFieldDecomposerPtrs_[patchi])
91  {
92  bf.set
93  (
94  patchi,
96  (
97  field.boundaryField()[boundaryAddressing_[patchi]],
98  procMesh_.boundary()[patchi],
99  resF.dimensionedInternalField(),
100  *patchFieldDecomposerPtrs_[patchi]
101  )
102  );
103  }
104  else if (isA<processorCyclicFvPatch>(procMesh_.boundary()[patchi]))
105  {
106  bf.set
107  (
108  patchi,
110  (
111  procMesh_.boundary()[patchi],
112  resF.dimensionedInternalField(),
114  (
115  field.internalField(),
116  *processorVolPatchFieldDecomposerPtrs_[patchi]
117  )
118  )
119  );
120  }
121  else if (isA<processorFvPatch>(procMesh_.boundary()[patchi]))
122  {
123  bf.set
124  (
125  patchi,
127  (
128  procMesh_.boundary()[patchi],
129  resF.dimensionedInternalField(),
131  (
132  field.internalField(),
133  *processorVolPatchFieldDecomposerPtrs_[patchi]
134  )
135  )
136  );
137  }
138  else if (allowUnknownPatchFields)
139  {
140  bf.set
141  (
142  patchi,
144  (
145  procMesh_.boundary()[patchi],
146  resF.dimensionedInternalField()
147  )
148  );
149  }
150  else
151  {
152  FatalErrorIn("fvFieldDecomposer::decomposeField()")
153  << "Unknown type." << abort(FatalError);
154  }
155  }
156 
157  // Create the field for the processor
158  return tresF;
159 }
160 
161 
162 template<class Type>
165 (
167 ) const
168 {
169  labelList mapAddr
170  (
172  (
173  faceAddressing_,
174  procMesh_.nInternalFaces()
175  )
176  );
177  forAll(mapAddr, i)
178  {
179  mapAddr[i] -= 1;
180  }
181 
182  // Create and map the internal field values
184  (
185  field.internalField(),
186  mapAddr
187  );
188 
189  // Problem with addressing when a processor patch picks up both internal
190  // faces and faces from cyclic boundaries. This is a bit of a hack, but
191  // I cannot find a better solution without making the internal storage
192  // mechanism for surfaceFields correspond to the one of faces in polyMesh
193  // (i.e. using slices)
194  Field<Type> allFaceField(field.mesh().nFaces());
195 
196  forAll(field.internalField(), i)
197  {
198  allFaceField[i] = field.internalField()[i];
199  }
200 
201  forAll(field.boundaryField(), patchi)
202  {
203  const Field<Type> & p = field.boundaryField()[patchi];
204 
205  const label patchStart = field.mesh().boundaryMesh()[patchi].start();
206 
207  forAll(p, i)
208  {
209  allFaceField[patchStart + i] = p[i];
210  }
211  }
212 
213 
214  // 1. Create the complete field with dummy patch fields
215  PtrList<fvsPatchField<Type> > patchFields(boundaryAddressing_.size());
216 
217  forAll(boundaryAddressing_, patchi)
218  {
219  patchFields.set
220  (
221  patchi,
223  (
225  procMesh_.boundary()[patchi],
227  )
228  );
229  }
230 
232  (
234  (
235  IOobject
236  (
237  field.name(),
238  procMesh_.time().timeName(),
239  procMesh_,
240  IOobject::NO_READ,
241  IOobject::NO_WRITE
242  ),
243  procMesh_,
244  field.dimensions(),
245  Field<Type>(field.internalField(), mapAddr),
246  patchFields
247  )
248  );
250 
251 
252  // 2. Change the fvsPatchFields to the correct type using a mapper
253  // constructor (with reference to the now correct internal field)
254 
256  GeometricBoundaryField& bf = resF.boundaryField();
257 
258  forAll(boundaryAddressing_, patchi)
259  {
260  if (patchFieldDecomposerPtrs_[patchi])
261  {
262  bf.set
263  (
264  patchi,
266  (
267  field.boundaryField()[boundaryAddressing_[patchi]],
268  procMesh_.boundary()[patchi],
269  resF.dimensionedInternalField(),
270  *patchFieldDecomposerPtrs_[patchi]
271  )
272  );
273  }
274  else if (isA<processorCyclicFvPatch>(procMesh_.boundary()[patchi]))
275  {
276  bf.set
277  (
278  patchi,
280  (
281  procMesh_.boundary()[patchi],
282  resF.dimensionedInternalField(),
284  (
285  allFaceField,
286  *processorSurfacePatchFieldDecomposerPtrs_[patchi]
287  )
288  )
289  );
290  }
291  else if (isA<processorFvPatch>(procMesh_.boundary()[patchi]))
292  {
293  bf.set
294  (
295  patchi,
297  (
298  procMesh_.boundary()[patchi],
299  resF.dimensionedInternalField(),
301  (
302  allFaceField,
303  *processorSurfacePatchFieldDecomposerPtrs_[patchi]
304  )
305  )
306  );
307  }
308  else
309  {
310  FatalErrorIn("fvFieldDecomposer::decomposeField()")
311  << "Unknown type." << abort(FatalError);
312  }
313  }
314 
315  // Create the field for the processor
316  return tresF;
317 }
318 
319 
320 template<class GeoField>
322 (
323  const PtrList<GeoField>& fields
324 ) const
325 {
326  forAll(fields, fieldI)
327  {
328  decomposeField(fields[fieldI])().write();
329  }
330 }
331 
332 
333 // ************************************************************************* //
This boundary condition provides an &#39;empty&#39; condition for reduced dimensions cases, i.e. 1- and 2-D geometries. Apply this condition to patches whose normal is aligned to geometric directions that do not constitue solution directions.
This boundary condition is not designed to be evaluated; it is assmued that the value is assigned via...
bool set(const label) const
Is element set.
Definition: PtrListI.H:107
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
This boundary condition enables processor communication across patches.
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
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
void decomposeFields(const PtrList< GeoField > &fields) const
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
InternalField & internalField()
Return internal field.
tmp< GeometricField< Type, fvPatchField, volMesh > > decomposeField(const GeometricField< Type, fvPatchField, volMesh > &field, const bool allowUnknownPatchFields=false) const
Decompose volume field.
const Mesh & mesh() const
Return mesh.
This boundary condition enables processor communication across cyclic patches.
runTime write()
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
volScalarField & p
Definition: createFields.H:51
#define forAll(list, i)
Definition: UList.H:421
label patchi
const dimensionSet & dimensions() const
Return dimensions.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Pre-declare SubField and related Field type.
Definition: Field.H:57
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const word & name() const
Return name.
Definition: IOobject.H:260
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Generic GeometricField class.
error FatalError
A List obtained as a section of another List.
Definition: SubList.H:53
Foam::calculatedFvsPatchField.
Foam::processorFvsPatchField.
A class for managing temporary objects.
Definition: PtrList.H:118
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
conserve internalField()+
Foam::processorCyclicFvsPatchField.