Loading...
Searching...
No Matches
setLayerPairing.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) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2020-2022 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
27Description
28 Remove a layer of cells and prepare addressing data
29
30\*---------------------------------------------------------------------------*/
31
33#include "polyMesh.H"
34#include "primitiveMesh.H"
35#include "polyTopoChange.H"
36#include "oppositeFace.H"
37#include "polyTopoChanger.H"
38
39// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
40
41bool Foam::layerAdditionRemoval::setLayerPairing() const
42{
43 // Note:
44 // This is also the most complex part of the topological change.
45 // Therefore it will be calculated here and stored as temporary
46 // data until the actual topological change, after which it will
47 // be cleared.
48
49 // Algorithm for point collapse
50 // 1) Go through the master cell layer and for every face of
51 // the face zone find the opposite face in the master cell.
52 // Check the direction of the opposite face and adjust as
53 // necessary. Check other faces to find an edge defining
54 // relative orientation of the two faces and adjust the face
55 // as necessary. Once the face is adjusted, record the
56 // addressing between the master and slave vertex layer.
57
58 const polyMesh& mesh = topoChanger().mesh();
59
60 const labelList& mc =
61 mesh.faceZones()[faceZoneID_.index()].masterCells();
62
63 const labelList& mf = mesh.faceZones()[faceZoneID_.index()];
64
65 const boolList& mfFlip =
66 mesh.faceZones()[faceZoneID_.index()].flipMap();
67
68 const faceList& faces = mesh.faces();
69 const cellList& cells = mesh.cells();
70
71 // Grab the local faces from the master zone
72 const faceList& mlf =
73 mesh.faceZones()[faceZoneID_.index()]().localFaces();
74
75 const labelList& meshPoints =
76 mesh.faceZones()[faceZoneID_.index()]().meshPoints();
77
78 // Create a list of points to collapse for every point of
79 // the master patch
80 if (pointsPairingPtr_ || facesPairingPtr_)
81 {
83 << "Problem with layer pairing data"
84 << abort(FatalError);
85 }
86
87 pointsPairingPtr_.reset(new labelList(meshPoints.size(), -1));
88 facesPairingPtr_.reset(new labelList(mf.size(), -1));
89
90 auto& ptc = *pointsPairingPtr_;
91 auto& ftc = *facesPairingPtr_;
92
93 if (debug > 1)
94 {
95 Pout<< "meshPoints: " << meshPoints << nl
96 << "localPoints: "
97 << mesh.faceZones()[faceZoneID_.index()]().localPoints()
98 << endl;
99 }
100
101 // For all faces, create the mapping
102 label nPointErrors = 0;
103 label nFaceErrors = 0;
104
105 forAll(mf, facei)
106 {
107 // Get the local master face
108 face curLocalFace = mlf[facei];
109
110 // Flip face based on flip index to recover original orientation
111 if (mfFlip[facei])
112 {
113 curLocalFace.flip();
114 }
115
116 // Get the opposing face from the master cell
117 oppositeFace lidFace =
118 cells[mc[facei]].opposingFace(mf[facei], faces);
119
120 if (!lidFace.good())
121 {
122 // This is not a valid layer; cannot continue
123 nFaceErrors++;
124 continue;
125 }
126
127 if (debug > 1)
128 {
129 Pout<< "curMasterFace: " << faces[mf[facei]] << nl
130 << "cell shape: " << mesh.cellShapes()[mc[facei]] << nl
131 << "curLocalFace: " << curLocalFace << nl
132 << "lidFace: " << lidFace
133 << " master index: " << lidFace.masterIndex()
134 << " oppositeIndex: " << lidFace.oppositeIndex() << endl;
135 }
136
137 // Grab the opposite face for face collapse addressing
138 ftc[facei] = lidFace.oppositeIndex();
139
140 // Using the local face insert the points into the lid list
141 forAll(curLocalFace, pointi)
142 {
143 const label clp = curLocalFace[pointi];
144
145 if (ptc[clp] == -1)
146 {
147 // Point not mapped yet. Insert the label
148 ptc[clp] = lidFace[pointi];
149 }
150 else
151 {
152 // Point mapped from some other face. Check the label
153 if (ptc[clp] != lidFace[pointi])
154 {
155 nPointErrors++;
156
157 if (debug > 1)
158 {
159 Pout<< "Topological error in cell layer pairing. "
160 << "This mesh is either topologically incorrect "
161 << "or the master face layer is not defined "
162 << "consistently. Please check the "
163 << "face zone flip map." << nl
164 << "First index: " << ptc[clp]
165 << " new index: " << lidFace[pointi] << endl;
166 }
167 }
168 }
169 }
170 }
171
172 if (returnReduceOr(nPointErrors || nFaceErrors))
173 {
174 clearAddressing();
175
176 return false;
177 }
178
179 // Valid layer
180 return true;
181}
182
183
184const Foam::labelList& Foam::layerAdditionRemoval::pointsPairing() const
185{
186 if (!pointsPairingPtr_)
187 {
189 << "Problem with layer pairing data for object " << name()
190 << abort(FatalError);
191 }
192
193 return *pointsPairingPtr_;
194}
195
196
197const Foam::labelList& Foam::layerAdditionRemoval::facesPairing() const
198{
199 if (!facesPairingPtr_)
200 {
202 << "Problem with layer pairing data for object " << name()
203 << abort(FatalError);
204 }
205
206 return *facesPairingPtr_;
207}
208
209
210// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
211
213(
214 pointField& motionPoints
215) const
216{
217 if (debug)
218 {
219 Pout<< "void layerAdditionRemoval::modifyMotionPoints("
220 << "pointField& motionPoints) const for object "
221 << name() << " : ";
222 }
223
224 if (debug)
225 {
226 Pout<< "No motion point adjustment" << endl;
227 }
228}
229
230
231// ************************************************************************* //
virtual void modifyMotionPoints(pointField &motionPoints) const
Modify motion points to comply with the topological change.
const polyTopoChanger & topoChanger() const
Return reference to morph engine.
const polyMesh & mesh() const
Return the mesh reference.
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition error.H:600
auto & name
const cellShapeList & cells
Namespace for handling debugging switches.
Definition debug.C:45
bool returnReduceOr(const bool value, const int communicator=UPstream::worldComm)
Perform logical (or) MPI Allreduce on a copy. Uses UPstream::reduceOr.
List< label > labelList
A List of labels.
Definition List.H:62
List< face > faceList
List of faces.
Definition faceListFwd.H:41
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition Ostream.H:519
List< cell > cellList
List of cell.
Definition cellListFwd.H:41
errorManip< error > abort(error &err)
Definition errorManip.H:139
List< bool > boolList
A List of bools.
Definition List.H:60
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
vectorField pointField
pointField is a vectorField.
constexpr char nl
The newline '\n' character (0x0a).
Definition Ostream.H:50
#define forAll(list, i)
Loop across all elements in list.
Definition stdFoam.H:299