sampledSurfacesTemplates.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-2021 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 "sampledSurfaces.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 #include "ListListOps.H"
30 #include "stringListOps.H"
31 
32 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33 
34 template<class Type>
36 Foam::functionObjects::sampledSurfaces::sampleLocalType
37 (
38  const label surfi,
39  const wordList& fieldNames,
40  HashPtrTable<interpolation<Type>>& interpolations
41 )
42 {
43  PtrList<Field<Type>> fieldTypeValues(fieldNames.size());
44 
45  const sampledSurface& s = operator[](surfi);
46 
47  forAll(fieldNames, fieldi)
48  {
49  const word& name = fieldNames[fieldi];
50 
52  {
53  const VolField<Type>& vf =
55 
56  if (s.interpolate())
57  {
58  if (!interpolations.found(name))
59  {
60  interpolations.insert
61  (
62  name,
64  (
65  interpolationScheme_,
66  vf
67  ).ptr()
68  );
69  }
70 
71  fieldTypeValues.set
72  (
73  fieldi,
74  s.interpolate(*interpolations[name]).ptr()
75  );
76  }
77  else
78  {
79  fieldTypeValues.set(fieldi, s.sample(vf).ptr());
80  }
81  }
83  {
84  const SurfaceField<Type>& sf =
86 
87  fieldTypeValues.set(fieldi, s.sample(sf).ptr());
88  }
89  }
90 
91  return fieldTypeValues;
92 }
93 
94 
95 template<class Type>
97 Foam::functionObjects::sampledSurfaces::sampleType
98 (
99  const label surfi,
100  const wordList& fieldNames,
101  HashPtrTable<interpolation<Type>>& interpolations
102 )
103 {
104  // Generate local samples
105  PtrList<Field<Type>> fieldTypeValues =
106  sampleLocalType<Type>(surfi, fieldNames, interpolations);
107 
108  if (Pstream::parRun())
109  {
110  // Collect values from all processors
111  PtrList<List<Field<Type>>> gatheredTypeValues(fieldNames.size());
112  forAll(fieldNames, fieldi)
113  {
114  if (fieldTypeValues.set(fieldi))
115  {
116  gatheredTypeValues.set
117  (
118  fieldi,
120  );
121  gatheredTypeValues[fieldi][Pstream::myProcNo()] =
122  fieldTypeValues[fieldi];
123  Pstream::gatherList(gatheredTypeValues[fieldi]);
124  }
125  }
126 
127  // Clear the local field values
128  fieldTypeValues.clear();
129  fieldTypeValues.resize(fieldNames.size());
130 
131  // Combine on the master
132  if (Pstream::master())
133  {
134  // Combine values into single field
135  forAll(fieldNames, fieldi)
136  {
137  if (gatheredTypeValues.set(fieldi))
138  {
139  fieldTypeValues.set
140  (
141  fieldi,
142  new Field<Type>
143  (
145  (
146  gatheredTypeValues[fieldi],
148  )
149  )
150  );
151  }
152  }
153 
154  // Renumber point data to correspond to merged points
155  forAll(fieldNames, fieldi)
156  {
157  if (fieldTypeValues.set(fieldi))
158  {
159  if
160  (
161  mergeList_[surfi].pointsMap.size()
162  == fieldTypeValues[fieldi].size()
163  )
164  {
165  Field<Type> f(fieldTypeValues[fieldi]);
166 
167  inplaceReorder(mergeList_[surfi].pointsMap, f);
168  f.setSize(mergeList_[surfi].points.size());
169 
170  fieldTypeValues.set(fieldi, new Field<Type>(f, true));
171  }
172  }
173  }
174  }
175  }
176 
177  return fieldTypeValues;
178 }
179 
180 
181 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Pre-declare SubField and related Field type.
Definition: Field.H:83
Generic GeometricField class.
A HashTable specialisation for hashing pointers.
Definition: HashPtrTable.H:67
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
bool set(const label) const
Is element set.
Definition: PtrListI.H:62
void resize(const label)
Alias for setSize(const label)
Definition: PtrListI.H:32
void clear()
Clear the PtrList, i.e. set size to zero deleting all the.
Definition: PtrList.C:198
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
const T & operator[](const label) const
Return element const reference.
Definition: UPtrListI.H:105
const word & name() const
Return the name of this functionObject.
const fvMesh & mesh_
Reference to the fvMesh.
Abstract base class for interpolation.
Definition: interpolation.H:55
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type and name.
bool foundObject(const word &name) const
Is the named Type in registry.
An abstract class for surfaces with sampling.
A class for handling words, derived from string.
Definition: word.H:62
volScalarField sf(fieldObject, mesh)
static List< word > fieldNames
Definition: globalFoam.H:46
const pointField & points
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
AccessType combine(const List< T > &, AccessOp aop=accessOp< T >())
Combines sublists into one big list.
Definition: ListListOps.C:34
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
void inplaceReorder(const labelUList &oldToNew, ListType &)
Inplace reorder the elements of a list.
labelList f(nPoints)
Operations on lists of strings.
Foam::surfaceFields.