faFieldDecomposerDecomposeFields.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2016-2017 Wikki Ltd
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "faFieldDecomposer.H"
29 #include "processorFaPatchField.H"
30 #include "processorFaePatchField.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
38 
39 template<class Type>
40 tmp<GeometricField<Type, faPatchField, areaMesh>>
42 (
43  const GeometricField<Type, faPatchField, areaMesh>& field
44 ) const
45 {
46  // Create and map the internal field values
47  Field<Type> internalField(field.internalField(), faceAddressing_);
48 
49  // Create and map the patch field values
50  PtrList<faPatchField<Type>> patchFields(boundaryAddressing_.size());
51 
52  forAll(boundaryAddressing_, patchi)
53  {
54  if (boundaryAddressing_[patchi] >= 0)
55  {
56  patchFields.set
57  (
58  patchi,
60  (
61  field.boundaryField()[boundaryAddressing_[patchi]],
62  procMesh_.boundary()[patchi],
64  *patchFieldDecomposerPtrs_[patchi]
65  )
66  );
67  }
68  else
69  {
70  patchFields.set
71  (
72  patchi,
73  new processorFaPatchField<Type>
74  (
75  procMesh_.boundary()[patchi],
77  Field<Type>
78  (
79  field.internalField(),
80  *processorAreaPatchFieldDecomposerPtrs_[patchi]
81  )
82  )
83  );
84  }
85  }
86 
87  // Create the field for the processor
88  return tmp<GeometricField<Type, faPatchField, areaMesh>>
89  (
90  new GeometricField<Type, faPatchField, areaMesh>
91  (
92  IOobject
93  (
94  field.name(),
95  procMesh_.time().timeName(),
96  procMesh_(),
99  ),
100  procMesh_,
101  field.dimensions(),
102  internalField,
103  patchFields
104  )
105  );
106 }
107 
108 
109 template<class Type>
110 tmp<GeometricField<Type, faePatchField, edgeMesh>>
112 (
113  const GeometricField<Type, faePatchField, edgeMesh>& field
114 ) const
115 {
116  labelList mapAddr
117  (
119  (
120  edgeAddressing_,
121  procMesh_.nInternalEdges()
122  )
123  );
124  forAll(mapAddr, i)
125  {
126  mapAddr[i] -= 1;
127  }
128 
129  // Create and map the internal field values
130  Field<Type> internalField
131  (
132  field.internalField(),
133  mapAddr
134  );
135 
136  // Problem with addressing when a processor patch picks up both internal
137  // edges and edges from cyclic boundaries. This is a bit of a hack, but
138  // I cannot find a better solution without making the internal storage
139  // mechanism for edgeFields correspond to the one of edges in polyMesh
140  // (i.e. using slices)
141  Field<Type> allEdgeField(field.mesh().nEdges());
142 
143  forAll(field.internalField(), i)
144  {
145  allEdgeField[i] = field.internalField()[i];
146  }
147 
148  forAll(field.boundaryField(), patchi)
149  {
150  const Field<Type> & p = field.boundaryField()[patchi];
151 
152  const label patchStart = field.mesh().boundary()[patchi].start();
153 
154  forAll(p, i)
155  {
156  allEdgeField[patchStart + i] = p[i];
157  }
158  }
159 
160  // Create and map the patch field values
161  PtrList<faePatchField<Type>> patchFields(boundaryAddressing_.size());
162 
163  forAll(boundaryAddressing_, patchi)
164  {
165  if (boundaryAddressing_[patchi] >= 0)
166  {
167  patchFields.set
168  (
169  patchi,
171  (
172  field.boundaryField()[boundaryAddressing_[patchi]],
173  procMesh_.boundary()[patchi],
175  *patchFieldDecomposerPtrs_[patchi]
176  )
177  );
178  }
179  else
180  {
181  patchFields.set
182  (
183  patchi,
184  new processorFaePatchField<Type>
185  (
186  procMesh_.boundary()[patchi],
188  Field<Type>
189  (
190  allEdgeField,
191  *processorEdgePatchFieldDecomposerPtrs_[patchi]
192  )
193  )
194  );
195  }
196  }
197 
198  // Create the field for the processor
199  return tmp<GeometricField<Type, faePatchField, edgeMesh>>
200  (
201  new GeometricField<Type, faePatchField, edgeMesh>
202  (
203  IOobject
204  (
205  field.name(),
206  procMesh_.time().timeName(),
207  procMesh_(),
210  ),
211  procMesh_,
212  field.dimensions(),
213  internalField,
214  patchFields
215  )
216  );
217 }
218 
219 
220 template<class GeoField>
222 (
223  const PtrList<GeoField>& fields
224 ) const
225 {
226  forAll(fields, fieldI)
227  {
228  decomposeField(fields[fieldI])().write();
229  }
230 }
231 
232 
233 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
234 
235 } // End namespace Foam
236 
237 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:130
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:71
p
volScalarField & p
Definition: createFieldRefs.H:8
processorFaPatchField.H
processorFaePatchField.H
Foam::faFieldDecomposer::decomposeField
tmp< GeometricField< Type, faPatchField, areaMesh > > decomposeField(const GeometricField< Type, faPatchField, areaMesh > &field) const
Decompose area field.
Foam::faPatchField::New
static tmp< faPatchField< Type > > New(const word &patchFieldType, const word &actualPatchType, const faPatch &, const DimensionedField< Type, areaMesh > &)
Definition: faPatchFieldNew.C:33
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::faePatchField::New
static tmp< faePatchField< Type > > New(const word &, const faPatch &, const DimensionedField< Type, edgeMesh > &)
Return a pointer to a new patchField created on freestore given.
Definition: faePatchFieldNew.C:33
Foam::DimensionedField::null
static const DimensionedField< Type, GeoMesh > & null()
Return a NullObjectRef DimensionedField.
Definition: DimensionedFieldI.H:33
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::faMesh::time
const Time & time() const
Return reference to time.
Definition: faMesh.C:897
field
rDeltaTY field()
Foam::List< label >::subList
SubList< label > subList
Declare type of subList.
Definition: List.H:114
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::faMesh::nInternalEdges
label nInternalEdges() const
Definition: faMesh.H:361
Foam::faMesh::boundary
const faBoundaryMesh & boundary() const
Return constant reference to boundary mesh.
Definition: faMesh.C:1019
Foam::vtk::write
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Definition: foamVtkOutputTemplates.C:35
Foam::faFieldDecomposer::decomposeFields
void decomposeFields(const PtrList< GeoField > &fields) const
Foam::IOobject::NO_READ
Definition: IOobject.H:123
fields
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:97
faFieldDecomposer.H