FreeFOAM The Cross-Platform CFD Toolkit
perfectInterface.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) 1991-2010 OpenCFD Ltd.
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 Description
25  Best thing is probably to look at attachDetach which does almost exactly
26  the same but for the geometric matching of points and face centres.
27 
28 \*---------------------------------------------------------------------------*/
29 
30 #include "perfectInterface.H"
32 #include <OpenFOAM/polyMesh.H>
35 #include <OpenFOAM/mapPolyMesh.H>
36 #include <OpenFOAM/matchPoints.H>
40 
41 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45  defineTypeNameAndDebug(perfectInterface, 0);
47  (
48  polyMeshModifier,
49  perfectInterface,
50  dictionary
51  );
52 }
53 
54 
55 // Tolerance used as fraction of minimum edge length.
56 const Foam::scalar Foam::perfectInterface::tol_ = 1E-3;
57 
58 
59 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
60 
61 Foam::pointField Foam::perfectInterface::calcFaceCentres
62 (
63  const primitivePatch& pp
64 )
65 {
66  const pointField& points = pp.points();
67 
68  pointField ctrs(pp.size());
69 
70  forAll(ctrs, patchFaceI)
71  {
72  ctrs[patchFaceI] = pp[patchFaceI].centre(points);
73  }
74 
75  return ctrs;
76 }
77 
78 
79 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
80 
81 // Construct from components
82 Foam::perfectInterface::perfectInterface
83 (
84  const word& name,
85  const label index,
86  const polyTopoChanger& mme,
87  const word& faceZoneName,
88  const word& masterPatchName,
89  const word& slavePatchName
90 )
91 :
92  polyMeshModifier(name, index, mme, true),
93  faceZoneID_(faceZoneName, mme.mesh().faceZones()),
94  masterPatchID_(masterPatchName, mme.mesh().boundaryMesh()),
95  slavePatchID_(slavePatchName, mme.mesh().boundaryMesh())
96 {}
97 
98 
99 // Construct from dictionary
100 Foam::perfectInterface::perfectInterface
101 (
102  const word& name,
103  const dictionary& dict,
104  const label index,
105  const polyTopoChanger& mme
106 )
107 :
108  polyMeshModifier(name, index, mme, readBool(dict.lookup("active"))),
109  faceZoneID_
110  (
111  dict.lookup("faceZoneName"),
112  mme.mesh().faceZones()
113  ),
114  masterPatchID_
115  (
116  dict.lookup("masterPatchName"),
117  mme.mesh().boundaryMesh()
118  ),
119  slavePatchID_
120  (
121  dict.lookup("slavePatchName"),
122  mme.mesh().boundaryMesh()
123  )
124 {}
125 
126 
127 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
128 
130 {}
131 
132 
133 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
134 
136 {
137  // If modifier is inactive, skip change
138  if (!active())
139  {
140  if (debug)
141  {
142  Pout<< "bool perfectInterface::changeTopology() const "
143  << "for object " << name() << " : "
144  << "Inactive" << endl;
145  }
146 
147  return false;
148  }
149  else
150  {
151  // I want topo change every time step.
152  return true;
153  }
154 }
155 
156 
158 {
159  if (debug)
160  {
161  Pout<< "bool perfectInterface::setRefinement(polyTopoChange&) const : "
162  << "for object " << name() << " : "
163  << "masterPatchID_:" << masterPatchID_
164  << " slavePatchID_:" << slavePatchID_
165  << " faceZoneID_:" << faceZoneID_ << endl;
166  }
167 
168  if
169  (
170  masterPatchID_.active()
171  && slavePatchID_.active()
172  && faceZoneID_.active()
173  )
174  {
175  const polyMesh& mesh = topoChanger().mesh();
176 
177  const polyBoundaryMesh& patches = mesh.boundaryMesh();
178 
179  const polyPatch& pp0 = patches[masterPatchID_.index()];
180  const polyPatch& pp1 = patches[slavePatchID_.index()];
181 
182  // Some aliases
183  const edgeList& edges0 = pp0.edges();
184  const pointField& pts0 = pp0.localPoints();
185  const pointField& pts1 = pp1.localPoints();
186  const labelList& meshPts0 = pp0.meshPoints();
187  const labelList& meshPts1 = pp1.meshPoints();
188 
189 
190  // Get local dimension as fraction of minimum edge length
191 
192  scalar minLen = GREAT;
193 
194  forAll(edges0, edgeI)
195  {
196  minLen = min(minLen, edges0[edgeI].mag(pts0));
197  }
198  scalar typDim = tol_*minLen;
199 
200  if (debug)
201  {
202  Pout<< "typDim:" << typDim << " edges0:" << edges0.size()
203  << " pts0:" << pts0.size() << " pts1:" << pts1.size()
204  << " pp0:" << pp0.size() << " pp1:" << pp1.size() << endl;
205  }
206 
207 
208  // Determine pointMapping in mesh point labels. Uses geometric
209  // comparison to find correspondence between patch points.
210 
211  labelList renumberPoints(mesh.points().size());
212  forAll(renumberPoints, i)
213  {
214  renumberPoints[i] = i;
215  }
216  {
217  labelList from1To0Points(pts1.size());
218 
219  bool matchOk = matchPoints
220  (
221  pts1,
222  pts0,
223  scalarField(pts1.size(), typDim), // tolerance
224  true, // verbose
225  from1To0Points
226  );
227 
228  if (!matchOk)
229  {
231  (
232  "perfectInterface::setRefinement(polyTopoChange& ref) const"
233  ) << "Points on patches " << pp0.name() << " and "
234  << pp1.name() << " do not match to within tolerance "
235  << typDim << exit(FatalError);
236  }
237 
238  forAll(pts1, i)
239  {
240  renumberPoints[meshPts1[i]] = meshPts0[from1To0Points[i]];
241  }
242  }
243 
244 
245 
246  // Calculate correspondence between patch faces
247 
248  labelList from0To1Faces(pp1.size());
249 
250  bool matchOk = matchPoints
251  (
252  calcFaceCentres(pp0),
253  calcFaceCentres(pp1),
254  scalarField(pp0.size(), typDim), // tolerance
255  true, // verbose
256  from0To1Faces
257  );
258 
259  if (!matchOk)
260  {
262  (
263  "perfectInterface::setRefinement(polyTopoChange& ref) const"
264  ) << "Face centres of patches " << pp0.name() << " and "
265  << pp1.name() << " do not match to within tolerance " << typDim
266  << exit(FatalError);
267  }
268 
269 
270 
271  // Now
272  // - renumber faces using pts1 (except patch1 faces)
273  // - remove patch1 faces. Remember cell label on owner side.
274  // - modify patch0 faces to be internal.
275 
276  // 1. Get faces to be renumbered
277  labelHashSet affectedFaces(2*pp1.size());
278  forAll(meshPts1, i)
279  {
280  label meshPointI = meshPts1[i];
281 
282  if (meshPointI != renumberPoints[meshPointI])
283  {
284  const labelList& pFaces = mesh.pointFaces()[meshPointI];
285 
286  forAll(pFaces, pFaceI)
287  {
288  affectedFaces.insert(pFaces[pFaceI]);
289  }
290  }
291  }
292  forAll(pp1, i)
293  {
294  affectedFaces.erase(pp1.start() + i);
295  }
296  // Remove patch0 from renumbered faces. Should not be nessecary since
297  // patch0 and 1 should not share any point (if created by mergeMeshing)
298  // so affectedFaces should not contain any patch0 faces but you can
299  // never be sure what the user is doing.
300  forAll(pp0, i)
301  {
302  if (affectedFaces.erase(pp0.start() + i))
303  {
304  WarningIn
305  (
306  "perfectInterface::setRefinement(polyTopoChange&) const"
307  ) << "Found face " << pp0.start() + i << " vertices "
308  << mesh.faces()[pp0.start() + i] << " whose points are"
309  << " used both by master patch " << pp0.name()
310  << " and slave patch " << pp1.name()
311  << endl;
312  }
313  }
314 
315 
316  // 2. Renumber (non patch0/1) faces.
317  for
318  (
319  labelHashSet::const_iterator iter = affectedFaces.begin();
320  iter != affectedFaces.end();
321  ++iter
322  )
323  {
324  label faceI = iter.key();
325 
326  const face& f = mesh.faces()[faceI];
327 
328  face newFace(f.size());
329 
330  forAll(newFace, fp)
331  {
332  newFace[fp] = renumberPoints[f[fp]];
333  }
334 
335  label nbr = -1;
336 
337  label patchI = -1;
338 
339  if (mesh.isInternalFace(faceI))
340  {
341  nbr = mesh.faceNeighbour()[faceI];
342  }
343  else
344  {
345  patchI = patches.whichPatch(faceI);
346  }
347 
348  label zoneID = mesh.faceZones().whichZone(faceI);
349 
350  bool zoneFlip = false;
351 
352  if (zoneID >= 0)
353  {
354  const faceZone& fZone = mesh.faceZones()[zoneID];
355 
356  zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)];
357  }
358 
359  ref.setAction
360  (
362  (
363  newFace, // modified face
364  faceI, // label of face being modified
365  mesh.faceOwner()[faceI], // owner
366  nbr, // neighbour
367  false, // face flip
368  patchI, // patch for face
369  false, // remove from zone
370  zoneID, // zone for face
371  zoneFlip // face flip in zone
372  )
373  );
374  }
375 
376 
377  // 3. Remove patch1 points
378  forAll(meshPts1, i)
379  {
380  label meshPointI = meshPts1[i];
381 
382  if (meshPointI != renumberPoints[meshPointI])
383  {
384  ref.setAction(polyRemovePoint(meshPointI));
385  }
386  }
387 
388 
389  // 4. Remove patch1 faces
390  forAll(pp1, i)
391  {
392  ref.setAction(polyRemoveFace(pp1.start() + i));
393  }
394 
395 
396  // 5. Modify patch0 faces for new points (not really nessecary; see
397  // comment above about patch1 and patch0 never sharing points) and
398  // becoming internal.
399  const boolList& mfFlip =
400  mesh.faceZones()[faceZoneID_.index()].flipMap();
401 
402  forAll(pp0, i)
403  {
404  label faceI = pp0.start() + i;
405 
406  const face& f = mesh.faces()[faceI];
407 
408  face newFace(f.size());
409 
410  forAll(newFace, fp)
411  {
412  newFace[fp] = renumberPoints[f[fp]];
413  }
414 
415  label own = mesh.faceOwner()[faceI];
416 
417  label pp1FaceI = pp1.start() + from0To1Faces[i];
418 
419  label nbr = mesh.faceOwner()[pp1FaceI];
420 
421  if (own < nbr)
422  {
423  ref.setAction
424  (
426  (
427  newFace, // modified face
428  faceI, // label of face being modified
429  own, // owner
430  nbr, // neighbour
431  false, // face flip
432  -1, // patch for face
433  false, // remove from zone
434  faceZoneID_.index(), // zone for face
435  mfFlip[i] // face flip in zone
436  )
437  );
438  }
439  else
440  {
441  ref.setAction
442  (
444  (
445  newFace.reverseFace(), // modified face
446  faceI, // label of face being modified
447  nbr, // owner
448  own, // neighbour
449  false, // face flip
450  -1, // patch for face
451  false, // remove from zone
452  faceZoneID_.index(), // zone for face
453  !mfFlip[i] // face flip in zone
454  )
455  );
456  }
457  }
458  }
459 }
460 
461 
463 {
464  // Update only my points. Nothing to be done here as points already
465  // shared by now.
466 }
467 
468 
470 {
471  // Mesh has changed topologically. Update local topological data
472  const polyMesh& mesh = topoChanger().mesh();
473 
474  faceZoneID_.update(mesh.faceZones());
475  masterPatchID_.update(mesh.boundaryMesh());
476  slavePatchID_.update(mesh.boundaryMesh());
477 }
478 
479 
481 {
482  os << nl << type() << nl
483  << name()<< nl
484  << faceZoneID_.name() << nl
485  << masterPatchID_.name() << nl
486  << slavePatchID_.name() << endl;
487 }
488 
489 
491 {
492  os << nl << name() << nl << token::BEGIN_BLOCK << nl
493 
494  << " type " << type()
496 
497  << " active " << active()
499 
500  << " faceZoneName " << faceZoneID_.name()
502 
503  << " masterPatchName " << masterPatchID_.name()
505 
506  << " slavePatchName " << slavePatchID_.name()
508 
509  << token::END_BLOCK << endl;
510 }
511 
512 
513 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
514 
515 
516 // * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
517 
518 
519 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
520 
521 
522 // ************************ vim: set sw=4 sts=4 et: ************************ //