Loading...
Searching...
No Matches
fvMeshAdder.H
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) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Class
28 Foam::fvMeshAdder
29
30Description
31 Adds two fvMeshes without using any polyMesh morphing.
32 Uses polyMeshAdder.
33
34SourceFiles
35 fvMeshAdder.C
36 fvMeshAdderTemplates.C
37
38\*---------------------------------------------------------------------------*/
39
40#ifndef fvMeshAdder_H
41#define fvMeshAdder_H
42
43#include "polyMeshAdder.H"
44#include "volFieldsFwd.H"
45#include "fvPatchFieldsFwd.H"
46#include "fvsPatchFieldsFwd.H"
47#include "fvPatchFieldMapper.H"
48#include "DimensionedField.H"
49
50// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51
52namespace Foam
53{
54
55// Forward Declarations
56class IOobject;
57class faceCoupleInfo;
58class IOobjectList;
59class fvMesh;
60class volMesh;
61class surfaceMesh;
64/*---------------------------------------------------------------------------*\
65 Class fvMeshAdder Declaration
66\*---------------------------------------------------------------------------*/
67
68class fvMeshAdder
69:
70 public polyMeshAdder
71{
72private:
73
74 // Private Member Functions
75
76 //- Calculate map from new patch faces to old patch faces. -1 where
77 // could not map.
78 static labelList calcPatchMap
79 (
80 const label oldStart,
81 const label oldSize,
82 const labelList& oldToNew,
83 const polyPatch& newPatch,
84 const label unmappedIndex
85 );
86
87 //- Update single volField.
88 template<class Type>
89 static void MapVolField
90 (
91 const mapAddedPolyMesh& meshMap,
92
95 const bool fullyMapped
96 );
97
98 //- Update single surfaceField.
99 template<class Type>
100 static void MapSurfaceField
101 (
102 const mapAddedPolyMesh& meshMap,
103
106 const bool fullyMapped
107 );
108
109 //- Update single dimensionedField.
110 template<class Type>
111 static void MapDimField
112 (
113 const mapAddedPolyMesh& meshMap,
114
116 const DimensionedField<Type, volMesh>& fldToAdd
117 );
118
119public:
120
121 // Declare name of the class and its debug switch
122 ClassName("fvMeshAdder");
123
124
125 // Member Functions
127 //- Inplace add mesh to fvMesh. Maps all stored fields. Returns map.
129 (
130 fvMesh& mesh0,
131 const fvMesh& mesh1,
132 const faceCoupleInfo& coupleInfo,
133 const bool validBoundary = true,
134 const bool fullyMapped = false
135 );
136
137 //- Map all volFields of Type.
138 // Optionally override mapping detection
139 // of unmapped values (e.g. used in fvMeshDistribute since it fixes
140 // up mapping itself)
141 template<class Type>
142 static void MapVolFields
143 (
144 const mapAddedPolyMesh&,
145 const fvMesh& mesh,
146 const fvMesh& meshToAdd,
147 const bool fullyMapped = false
148 );
149
150 //- Map all surfaceFields of Type
151 template<class Type>
152 static void MapSurfaceFields
153 (
154 const mapAddedPolyMesh&,
155 const fvMesh& mesh,
156 const fvMesh& meshToAdd,
157 const bool fullyMapped = false
158 );
159
160 //- Map all DimensionedFields of Type
161 template<class Type>
162 static void MapDimFields
163 (
164 const mapAddedPolyMesh&,
165 const fvMesh& mesh,
166 const fvMesh& meshToAdd
167 );
168
169
170 // Multi-mesh merging. Expects same patches (apart from processor
171 // patches). Adds to element 0 of flds, meshes
172
173 //- In-place add to fvMeshes[myProci]. Stitch boundary faces.
175 (
176 const label myProci, // index of mesh to modify (== mesh_)
177 UPtrList<fvMesh>& fvMeshes,
178 const labelList& oldFaceOwner, // face owner for myProci mesh
179
180 // Coupling info
181 const labelListList& localBoundaryFace,
182 const labelListList& remoteFaceProc,
183 const labelListList& remoteBoundaryFace,
184
185 labelListList& constructPatchMap,
186 labelListList& constructCellMap,
187 labelListList& constructFaceMap,
188 labelListList& constructPointMap
189 );
190
191 //- Update single dimensionedField.
192 template<class Type>
193 static void MapDimField
194 (
196 const labelListList& cellProcAddressing,
197 const bool fullyMapped
198 );
199
200 //- Update single volField.
201 template<class Type>
202 static void MapVolField
203 (
205 const labelList& oldPatchStarts0,
206 const labelList& oldPatchSizes0,
207 const labelListList& patchProcAddressing,
208 const labelListList& cellProcAddressing,
209 const labelListList& faceProcAddressing,
210 const bool fullyMapped
211 );
212
213 //- Update single surfaceField.
214 template<class Type>
215 static void MapSurfaceField
216 (
218 const labelList& oldFaceOwner0,
219 const labelList& oldPatchStarts0,
220 const labelList& oldPatchSizes0,
221 const labelListList& patchProcAddressing,
222 const labelListList& cellProcAddressing,
223 const labelListList& faceProcAddressing,
224 const bool fullyMapped
225 );
226
227
228 //- Map all dimensionedField of Type.
229 // Optionally override mapping detection
230 // of unmapped values (e.g. used in fvMeshDistribute since it
231 // fixes up mapping itself)
232 template<class Type>
233 static void MapDimFields
234 (
236 const labelListList& cellProcAddressing,
237 const bool fullyMapped = false
238 );
239
240 //- Map all volFields of Type.
241 // Optionally override mapping detection
242 // of unmapped values (e.g. used in fvMeshDistribute since it
243 // fixes up mapping itself)
244 template<class Type>
245 static void MapVolFields
246 (
248 const labelList& oldPatchStarts0,
249 const labelList& oldPatchSizes0,
250 const labelListList& patchProcAddressing,
251 const labelListList& cellProcAddressing,
252 const labelListList& faceProcAddressing,
253 const labelListList& pointProcAddressing,
254 const bool fullyMapped = false
255 );
256
257 //- Map all surfaceFields of Type
258 template<class Type>
259 static void MapSurfaceFields
260 (
262 const labelList& oldFaceOwner0,
263 const labelList& oldPatchStarts0,
264 const labelList& oldPatchSizes0,
265 const labelListList& patchProcAddressing,
266 const labelListList& cellProcAddressing,
267 const labelListList& faceProcAddressing,
268 const labelListList& pointProcAddressing,
269 const bool fullyMapped = false
270 );
271};
272
273
274// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
275
276} // End namespace Foam
277
278// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
279
280#ifdef NoRepository
281 #include "fvMeshAdderTemplates.C"
282#endif
283
284// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
285
286#endif
287
288// ************************************************************************* //
Info<< nl;Info<< "Write faMesh in vtk format:"<< nl;{ vtk::uindirectPatchWriter writer(aMesh.patch(), fileName(aMesh.time().globalPath()/vtkBaseFileName));writer.writeGeometry();globalIndex procAddr(aMesh.nFaces());labelList cellIDs;if(UPstream::master()) { cellIDs.resize(procAddr.totalSize());for(const labelRange &range :procAddr.ranges()) { auto slice=cellIDs.slice(range);slice=identity(range);} } writer.beginCellData(4);writer.writeProcIDs();writer.write("cellID", cellIDs);writer.write("area", aMesh.S().field());writer.write("normal", aMesh.faceAreaNormals());writer.beginPointData(1);writer.write("normal", aMesh.pointAreaNormals());Info<< " "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.time().globalPath()/(vtkBaseFileName+"-edges")));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
List of IOobjects with searching and retrieving facilities. Implemented as a HashTable,...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition IOobject.H:191
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
Definition UPtrList.H:101
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition autoPtr.H:65
Container for information needed to couple to meshes. When constructed from two meshes and a geometri...
Adds two fvMeshes without using any polyMesh morphing. Uses polyMeshAdder.
Definition fvMeshAdder.H:66
static void MapSurfaceFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd, const bool fullyMapped=false)
Map all surfaceFields of Type.
static autoPtr< mapAddedPolyMesh > add(fvMesh &mesh0, const fvMesh &mesh1, const faceCoupleInfo &coupleInfo, const bool validBoundary=true, const bool fullyMapped=false)
Inplace add mesh to fvMesh. Maps all stored fields. Returns map.
Definition fvMeshAdder.C:67
static void MapDimFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all DimensionedFields of Type.
static void MapVolFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd, const bool fullyMapped=false)
Map all volFields of Type.
ClassName("fvMeshAdder")
Mesh data needed to do the Finite Volume discretisation.
Definition fvMesh.H:85
Class containing mesh-to-mesh mapping information after a mesh addition where we add a mesh ('added m...
Adds two meshes without using any polyMesh morphing.
A patch is a list of labels that address the faces in the global face list.
Definition polyPatch.H:73
Mesh data needed to do the Finite Volume discretisation.
Definition surfaceMesh.H:47
Mesh data needed to do the Finite Volume discretisation.
Definition volMesh.H:47
#define ClassName(TypeNameString)
Add typeName information from argument TypeNameString to a class.
Definition className.H:74
dynamicFvMesh & mesh
Foam::PtrList< Foam::fvMesh > meshes(regionNames.size())
Namespace for OpenFOAM.
List< labelList > labelListList
List of labelList.
Definition labelList.H:38
List< label > labelList
A List of labels.
Definition List.H:62
Forwards and collection of common volume field types.