createBaffles.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) 2011-2016 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 Application
25  createBaffles
26 
27 Description
28  Makes internal faces into boundary faces. Does not duplicate points, unlike
29  mergeOrSplitBaffles.
30 
31  Note: if any coupled patch face is selected for baffling the opposite
32  member has to be selected for baffling as well.
33 
34  - if the patch already exists will not override it nor its fields
35  - if the patch does not exist it will be created together with 'calculated'
36  patchfields unless the field is mentioned in the patchFields section.
37  - any 0-sized patches (since faces have been moved out) will get removed
38 
39 \*---------------------------------------------------------------------------*/
40 
41 #include "argList.H"
42 #include "Time.H"
43 #include "polyTopoChange.H"
44 #include "polyModifyFace.H"
45 #include "polyAddFace.H"
46 #include "ReadFields.H"
47 #include "volFields.H"
48 #include "surfaceFields.H"
49 #include "fvMeshMapper.H"
50 #include "faceSelection.H"
51 
52 #include "fvMeshTools.H"
53 
54 using namespace Foam;
55 
56 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
57 
58 label addPatch
59 (
60  fvMesh& mesh,
61  const word& patchName,
62  const word& groupName,
63  const dictionary& patchDict
64 )
65 {
66  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
67 
68  if (pbm.findPatchID(patchName) == -1)
69  {
70  autoPtr<polyPatch> ppPtr
71  (
73  (
74  patchName,
75  patchDict,
76  0,
77  pbm
78  )
79  );
80  polyPatch& pp = ppPtr();
81 
82  if (!groupName.empty() && !pp.inGroup(groupName))
83  {
84  pp.inGroups().append(groupName);
85  }
86 
87 
88  // Add patch, create calculated everywhere
90  (
91  mesh,
92  pp,
93  dictionary(), // do not set specialised patchFields
95  true // parallel sync'ed addition
96  );
97  }
98  else
99  {
100  Info<< "Patch '" << patchName
101  << "' already exists. Only "
102  << "moving patch faces - type will remain the same"
103  << endl;
104  }
105 
106  return pbm.findPatchID(patchName);
107 }
108 
109 
110 // Filter out the empty patches.
111 void filterPatches(fvMesh& mesh, const HashSet<word>& addedPatchNames)
112 {
113  // Remove any zero-sized ones. Assumes
114  // - processor patches are already only there if needed
115  // - all other patches are available on all processors
116  // - but coupled ones might still be needed, even if zero-size
117  // (e.g. processorCyclic)
118  // See also logic in createPatch.
119  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
120 
121  labelList oldToNew(pbm.size(), -1);
122  label newPatchi = 0;
123  forAll(pbm, patchi)
124  {
125  const polyPatch& pp = pbm[patchi];
126 
127  if (!isA<processorPolyPatch>(pp))
128  {
129  if
130  (
131  isA<coupledPolyPatch>(pp)
132  || returnReduce(pp.size(), sumOp<label>())
133  || addedPatchNames.found(pp.name())
134  )
135  {
136  // Coupled (and unknown size) or uncoupled and used
137  oldToNew[patchi] = newPatchi++;
138  }
139  }
140  }
141 
142  forAll(pbm, patchi)
143  {
144  const polyPatch& pp = pbm[patchi];
145 
146  if (isA<processorPolyPatch>(pp))
147  {
148  oldToNew[patchi] = newPatchi++;
149  }
150  }
151 
152 
153  const label nKeepPatches = newPatchi;
154 
155  // Shuffle unused ones to end
156  if (nKeepPatches != pbm.size())
157  {
158  Info<< endl
159  << "Removing zero-sized patches:" << endl << incrIndent;
160 
161  forAll(oldToNew, patchi)
162  {
163  if (oldToNew[patchi] == -1)
164  {
165  Info<< indent << pbm[patchi].name()
166  << " type " << pbm[patchi].type()
167  << " at position " << patchi << endl;
168  oldToNew[patchi] = newPatchi++;
169  }
170  }
171  Info<< decrIndent;
172 
173  fvMeshTools::reorderPatches(mesh, oldToNew, nKeepPatches, true);
174  Info<< endl;
175  }
176 }
177 
178 
179 void modifyOrAddFace
180 (
181  polyTopoChange& meshMod,
182  const face& f,
183  const label facei,
184  const label own,
185  const bool flipFaceFlux,
186  const label newPatchi,
187  const label zoneID,
188  const bool zoneFlip,
189 
190  PackedBoolList& modifiedFace
191 )
192 {
193  if (!modifiedFace[facei])
194  {
195  // First usage of face. Modify.
196  meshMod.setAction
197  (
199  (
200  f, // modified face
201  facei, // label of face
202  own, // owner
203  -1, // neighbour
204  flipFaceFlux, // face flip
205  newPatchi, // patch for face
206  false, // remove from zone
207  zoneID, // zone for face
208  zoneFlip // face flip in zone
209  )
210  );
211  modifiedFace[facei] = 1;
212  }
213  else
214  {
215  // Second or more usage of face. Add.
216  meshMod.setAction
217  (
219  (
220  f, // modified face
221  own, // owner
222  -1, // neighbour
223  -1, // master point
224  -1, // master edge
225  facei, // master face
226  flipFaceFlux, // face flip
227  newPatchi, // patch for face
228  zoneID, // zone for face
229  zoneFlip // face flip in zone
230  )
231  );
232  }
233 }
234 
235 
236 // Create faces for fZone faces. Usually newMasterPatches, newSlavePatches
237 // only size one but can be more for duplicate baffle sets
238 void createFaces
239 (
240  const bool internalFacesOnly,
241  const fvMesh& mesh,
242  const faceZone& fZone,
243  const labelList& newMasterPatches,
244  const labelList& newSlavePatches,
245  polyTopoChange& meshMod,
246  PackedBoolList& modifiedFace,
247  label& nModified
248 )
249 {
250  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
251 
252  forAll(newMasterPatches, i)
253  {
254  // Pass 1. Do selected side of zone
255  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
256 
257  for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
258  {
259  label zoneFacei = fZone.whichFace(facei);
260 
261  if (zoneFacei != -1)
262  {
263  if (!fZone.flipMap()[zoneFacei])
264  {
265  // Use owner side of face
266  modifyOrAddFace
267  (
268  meshMod,
269  mesh.faces()[facei], // modified face
270  facei, // label of face
271  mesh.faceOwner()[facei],// owner
272  false, // face flip
273  newMasterPatches[i], // patch for face
274  fZone.index(), // zone for face
275  false, // face flip in zone
276  modifiedFace // modify or add status
277  );
278  }
279  else
280  {
281  // Use neighbour side of face.
282  // To keep faceZone pointing out of original neighbour
283  // we don't need to set faceFlip since that cell
284  // now becomes the owner
285  modifyOrAddFace
286  (
287  meshMod,
288  mesh.faces()[facei].reverseFace(), // modified face
289  facei, // label of face
290  mesh.faceNeighbour()[facei],// owner
291  true, // face flip
292  newMasterPatches[i], // patch for face
293  fZone.index(), // zone for face
294  false, // face flip in zone
295  modifiedFace // modify or add status
296  );
297  }
298 
299  nModified++;
300  }
301  }
302 
303 
304  // Pass 2. Do other side of zone
305  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
306 
307  for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
308  {
309  label zoneFacei = fZone.whichFace(facei);
310 
311  if (zoneFacei != -1)
312  {
313  if (!fZone.flipMap()[zoneFacei])
314  {
315  // Use neighbour side of face
316  modifyOrAddFace
317  (
318  meshMod,
319  mesh.faces()[facei].reverseFace(), // modified face
320  facei, // label of face
321  mesh.faceNeighbour()[facei], // owner
322  true, // face flip
323  newSlavePatches[i], // patch for face
324  fZone.index(), // zone for face
325  true, // face flip in zone
326  modifiedFace // modify or add
327  );
328  }
329  else
330  {
331  // Use owner side of face
332  modifyOrAddFace
333  (
334  meshMod,
335  mesh.faces()[facei], // modified face
336  facei, // label of face
337  mesh.faceOwner()[facei],// owner
338  false, // face flip
339  newSlavePatches[i], // patch for face
340  fZone.index(), // zone for face
341  true, // face flip in zone
342  modifiedFace // modify or add status
343  );
344  }
345  }
346  }
347 
348 
349  // Modify any boundary faces
350  // ~~~~~~~~~~~~~~~~~~~~~~~~~
351 
352  // Normal boundary:
353  // - move to new patch. Might already be back-to-back baffle
354  // you want to add cyclic to. Do warn though.
355  //
356  // Processor boundary:
357  // - do not move to cyclic
358  // - add normal patches though.
359 
360  // For warning once per patch.
361  labelHashSet patchWarned;
362 
363  forAll(pbm, patchi)
364  {
365  const polyPatch& pp = pbm[patchi];
366 
367  label newPatchi = newMasterPatches[i];
368 
369  if (pp.coupled() && pbm[newPatchi].coupled())
370  {
371  // Do not allow coupled faces to be moved to different
372  // coupled patches.
373  }
374  else if (pp.coupled() || !internalFacesOnly)
375  {
376  forAll(pp, i)
377  {
378  label facei = pp.start()+i;
379 
380  label zoneFacei = fZone.whichFace(facei);
381 
382  if (zoneFacei != -1)
383  {
384  if (patchWarned.insert(patchi))
385  {
387  << "Found boundary face (in patch "
388  << pp.name()
389  << ") in faceZone " << fZone.name()
390  << " to convert to baffle patch "
391  << pbm[newPatchi].name()
392  << endl
393  << " Run with -internalFacesOnly option"
394  << " if you don't wish to convert"
395  << " boundary faces." << endl;
396  }
397 
398  modifyOrAddFace
399  (
400  meshMod,
401  mesh.faces()[facei], // modified face
402  facei, // label of face
403  mesh.faceOwner()[facei], // owner
404  false, // face flip
405  newPatchi, // patch for face
406  fZone.index(), // zone for face
407  fZone.flipMap()[zoneFacei], // face flip in zone
408  modifiedFace // modify or add
409  );
410  nModified++;
411  }
412  }
413  }
414  }
415  }
416 }
417 
418 
419 int main(int argc, char *argv[])
420 {
422  (
423  "Makes internal faces into boundary faces.\n"
424  "Does not duplicate points."
425  );
426  #include "addDictOption.H"
427  #include "addOverwriteOption.H"
428  #include "addDictOption.H"
429  #include "addRegionOption.H"
430  #include "setRootCase.H"
431  #include "createTime.H"
432  runTime.functionObjects().off();
433  #include "createNamedMesh.H"
434 
435 
436  const bool overwrite = args.optionFound("overwrite");
437 
438  const word oldInstance = mesh.pointsInstance();
439 
440  const word dictName("createBafflesDict");
441  #include "setSystemMeshDictionaryIO.H"
442 
443  Switch internalFacesOnly(false);
444 
445  Switch noFields(false);
446 
447  PtrList<faceSelection> selectors;
448  {
449  Info<< "Reading baffle criteria from " << dictName << nl << endl;
451 
452  dict.lookup("internalFacesOnly") >> internalFacesOnly;
453  noFields = dict.lookupOrDefault("noFields", false);
454 
455  const dictionary& selectionsDict = dict.subDict("baffles");
456 
457  label n = 0;
458  forAllConstIter(dictionary, selectionsDict, iter)
459  {
460  if (iter().isDict())
461  {
462  n++;
463  }
464  }
465  selectors.setSize(n);
466  n = 0;
467  forAllConstIter(dictionary, selectionsDict, iter)
468  {
469  if (iter().isDict())
470  {
471  selectors.set
472  (
473  n++,
474  faceSelection::New(iter().keyword(), mesh, iter().dict())
475  );
476  }
477  }
478  }
479 
480 
481  if (internalFacesOnly)
482  {
483  Info<< "Not converting faces on non-coupled patches." << nl << endl;
484  }
485 
486 
487  // Read objects in time directory
488  IOobjectList objects(mesh, runTime.timeName());
489 
490  // Read vol fields.
491  Info<< "Reading geometric fields" << nl << endl;
493  if (!noFields) ReadFields(mesh, objects, vsFlds);
494 
496  if (!noFields) ReadFields(mesh, objects, vvFlds);
497 
499  if (!noFields) ReadFields(mesh, objects, vstFlds);
500 
501  PtrList<volSymmTensorField> vsymtFlds;
502  if (!noFields) ReadFields(mesh, objects, vsymtFlds);
503 
505  if (!noFields) ReadFields(mesh, objects, vtFlds);
506 
507  // Read surface fields.
508 
510  if (!noFields) ReadFields(mesh, objects, ssFlds);
511 
513  if (!noFields) ReadFields(mesh, objects, svFlds);
514 
516  if (!noFields) ReadFields(mesh, objects, sstFlds);
517 
519  if (!noFields) ReadFields(mesh, objects, ssymtFlds);
520 
522  if (!noFields) ReadFields(mesh, objects, stFlds);
523 
524 
525 
526 
527  // Creating (if necessary) faceZones
528  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
529 
530  forAll(selectors, selectorI)
531  {
532  const word& name = selectors[selectorI].name();
533 
534  if (mesh.faceZones().findZoneID(name) == -1)
535  {
536  mesh.faceZones().clearAddressing();
537  label sz = mesh.faceZones().size();
538 
539  labelList addr(0);
540  boolList flip(0);
541  mesh.faceZones().setSize(sz+1);
542  mesh.faceZones().set
543  (
544  sz,
545  new faceZone(name, addr, flip, sz, mesh.faceZones())
546  );
547  }
548  }
549 
550 
551  // Select faces
552  // ~~~~~~~~~~~~
553 
554  //- Per face zoneID it is in and flip status.
555  labelList faceToZoneID(mesh.nFaces(), -1);
556  boolList faceToFlip(mesh.nFaces(), false);
557  forAll(selectors, selectorI)
558  {
559  const word& name = selectors[selectorI].name();
560  label zoneID = mesh.faceZones().findZoneID(name);
561 
562  selectors[selectorI].select(zoneID, faceToZoneID, faceToFlip);
563  }
564 
565  // Add faces to faceZones
566  labelList nFaces(mesh.faceZones().size(), 0);
567  forAll(faceToZoneID, facei)
568  {
569  label zoneID = faceToZoneID[facei];
570  if (zoneID != -1)
571  {
572  nFaces[zoneID]++;
573  }
574  }
575 
576  forAll(selectors, selectorI)
577  {
578  const word& name = selectors[selectorI].name();
579  label zoneID = mesh.faceZones().findZoneID(name);
580 
581  label& n = nFaces[zoneID];
582  labelList addr(n);
583  boolList flip(n);
584  n = 0;
585  forAll(faceToZoneID, facei)
586  {
587  label zone = faceToZoneID[facei];
588  if (zone == zoneID)
589  {
590  addr[n] = facei;
591  flip[n] = faceToFlip[facei];
592  n++;
593  }
594  }
595 
596  Info<< "Created zone " << name
597  << " at index " << zoneID
598  << " with " << n << " faces" << endl;
599 
600  mesh.faceZones().set
601  (
602  zoneID,
603  new faceZone(name, addr, flip, zoneID, mesh.faceZones())
604  );
605  }
606 
607 
608 
609  // Count patches to add
610  // ~~~~~~~~~~~~~~~~~~~~
611  HashSet<word> bafflePatches;
612  {
613  forAll(selectors, selectorI)
614  {
615  const dictionary& dict = selectors[selectorI].dict();
616 
617  if (dict.found("patches"))
618  {
619  const dictionary& patchSources = dict.subDict("patches");
620  forAllConstIter(dictionary, patchSources, iter)
621  {
622  const word patchName(iter().dict()["name"]);
623  bafflePatches.insert(patchName);
624  }
625  }
626  else
627  {
628  const word masterName = selectors[selectorI].name() + "_master";
629  bafflePatches.insert(masterName);
630  const word slaveName = selectors[selectorI].name() + "_slave";
631  bafflePatches.insert(slaveName);
632  }
633  }
634  }
635 
636 
637  // Create baffles
638  // ~~~~~~~~~~~~~~
639  // Is done in multiple steps
640  // - create patches with 'calculated' patchFields
641  // - move faces into these patches
642  // - change the patchFields to the wanted type
643  // This order is done so e.g. fixedJump works:
644  // - you cannot create patchfields at the same time as patches since
645  // they do an evaluate upon construction
646  // - you want to create the patchField only after you have faces
647  // so you don't get the 'create-from-nothing' mapping problem.
648 
649 
650  // Pass 1: add patches
651  // ~~~~~~~~~~~~~~~~~~~
652 
653  //HashSet<word> addedPatches;
654  {
655  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
656  forAll(selectors, selectorI)
657  {
658  const dictionary& dict = selectors[selectorI].dict();
659  const word& groupName = selectors[selectorI].name();
660 
661  if (dict.found("patches"))
662  {
663  const dictionary& patchSources = dict.subDict("patches");
664  forAllConstIter(dictionary, patchSources, iter)
665  {
666  const word patchName(iter().dict()["name"]);
667 
668  if (pbm.findPatchID(patchName) == -1)
669  {
670  dictionary patchDict = iter().dict();
671  patchDict.set("nFaces", 0);
672  patchDict.set("startFace", 0);
673 
674  // Note: do not set coupleGroup if constructed from
675  // baffles so you have freedom specifying it
676  // yourself.
677  //patchDict.set("coupleGroup", groupName);
678 
679  addPatch(mesh, patchName, groupName, patchDict);
680  }
681  else
682  {
683  Info<< "Patch '" << patchName
684  << "' already exists. Only "
685  << "moving patch faces - type will remain the same"
686  << endl;
687  }
688  }
689  }
690  else
691  {
692  const dictionary& patchSource = dict.subDict("patchPairs");
693  const word masterName = groupName + "_master";
694  const word slaveName = groupName + "_slave";
695 
696  word groupNameMaster = groupName;
697  word groupNameSlave = groupName;
698 
699 
700  dictionary patchDictMaster(patchSource);
701  patchDictMaster.set("nFaces", 0);
702  patchDictMaster.set("startFace", 0);
703  patchDictMaster.set("coupleGroup", groupName);
704 
705  dictionary patchDictSlave(patchDictMaster);
706 
707  // Note: This is added for the particular case where we want
708  // master and slave in different groupNames
709  // (ie 3D thermal baffles)
710 
711  Switch sameGroup
712  (
713  patchSource.lookupOrDefault("sameGroup", true)
714  );
715  if (!sameGroup)
716  {
717  groupNameMaster = groupName + "Group_master";
718  groupNameSlave = groupName + "Group_slave";
719  patchDictMaster.set("coupleGroup", groupNameMaster);
720  patchDictSlave.set("coupleGroup", groupNameSlave);
721  }
722 
723  addPatch(mesh, masterName, groupNameMaster, patchDictMaster);
724  addPatch(mesh, slaveName, groupNameSlave, patchDictSlave);
725  }
726  }
727  }
728 
729 
730  // Make sure patches and zoneFaces are synchronised across couples
731  mesh.boundaryMesh().checkParallelSync(true);
732  mesh.faceZones().checkParallelSync(true);
733 
734 
735 
736  // Mesh change container
737  polyTopoChange meshMod(mesh);
738 
739  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
740 
741 
742  // Do the actual changes. Note:
743  // - loop in incrementing face order (not necessary if faceZone ordered).
744  // Preserves any existing ordering on patch faces.
745  // - two passes, do non-flip faces first and flip faces second. This
746  // guarantees that when e.g. creating a cyclic all faces from one
747  // side come first and faces from the other side next.
748 
749  // Whether first use of face (modify) or consecutive (add)
750  PackedBoolList modifiedFace(mesh.nFaces());
751  label nModified = 0;
752 
753  forAll(selectors, selectorI)
754  {
755  const word& name = selectors[selectorI].name();
756  label zoneID = mesh.faceZones().findZoneID(name);
757  const faceZone& fZone = mesh.faceZones()[zoneID];
758 
759  const dictionary& dict = selectors[selectorI].dict();
760 
761  DynamicList<label> newMasterPatches;
762  DynamicList<label> newSlavePatches;
763 
764  if (dict.found("patches"))
765  {
766  const dictionary& patchSources = dict.subDict("patches");
767 
768  bool master = true;
769  forAllConstIter(dictionary, patchSources, iter)
770  {
771  const word patchName(iter().dict()["name"]);
772  label patchi = pbm.findPatchID(patchName);
773  if (master)
774  {
775  newMasterPatches.append(patchi);
776  }
777  else
778  {
779  newSlavePatches.append(patchi);
780  }
781  master = !master;
782  }
783  }
784  else
785  {
786  const word masterName = selectors[selectorI].name() + "_master";
787  newMasterPatches.append(pbm.findPatchID(masterName));
788 
789  const word slaveName = selectors[selectorI].name() + "_slave";
790  newSlavePatches.append(pbm.findPatchID(slaveName));
791  }
792 
793 
794 
795  createFaces
796  (
797  internalFacesOnly,
798  mesh,
799  fZone,
800  newMasterPatches,
801  newSlavePatches,
802  meshMod,
803  modifiedFace,
804  nModified
805  );
806  }
807 
808 
809  Info<< "Converted " << returnReduce(nModified, sumOp<label>())
810  << " faces into boundary faces in patches "
811  << bafflePatches.sortedToc() << nl << endl;
812 
813  if (!overwrite)
814  {
815  runTime++;
816  }
817 
818  // Change the mesh. Change points directly (no inflation).
819  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
820 
821  // Update fields
822  mesh.updateMesh(map);
823 
824 
825 
826  // Correct boundary faces mapped-out-of-nothing.
827  // This is just a hack to correct the value field.
828  {
829  fvMeshMapper mapper(mesh, map);
830  bool hasWarned = false;
831 
832  forAllConstIter(HashSet<word>, bafflePatches, iter)
833  {
834  label patchi = mesh.boundaryMesh().findPatchID(iter.key());
835 
836  const fvPatchMapper& pm = mapper.boundaryMap()[patchi];
837 
838  if (pm.sizeBeforeMapping() == 0)
839  {
840  if (!hasWarned)
841  {
842  hasWarned = true;
844  << "Setting field on boundary faces to zero." << endl
845  << "You might have to edit these fields." << endl;
846  }
847 
848  fvMeshTools::zeroPatchFields(mesh, patchi);
849  }
850  }
851  }
852 
853 
854  // Pass 2: change patchFields
855  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
856 
857  {
858  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
859  forAll(selectors, selectorI)
860  {
861  const dictionary& dict = selectors[selectorI].dict();
862  if (dict.found("patches"))
863  {
864  const dictionary& patchSources = dict.subDict("patches");
865 
866  forAllConstIter(dictionary, patchSources, iter)
867  {
868  const word patchName(iter().dict()["name"]);
869  label patchi = pbm.findPatchID(patchName);
870 
871  if (iter().dict().found("patchFields"))
872  {
873  const dictionary& patchFieldsDict =
874  iter().dict().subDict
875  (
876  "patchFields"
877  );
878 
879  fvMeshTools::setPatchFields
880  (
881  mesh,
882  patchi,
883  patchFieldsDict
884  );
885  }
886  }
887  }
888  else
889  {
890  const dictionary& patchSource = dict.subDict("patchPairs");
891 
892  Switch sameGroup
893  (
894  patchSource.lookupOrDefault("sameGroup", true)
895  );
896 
897  const word& groupName = selectors[selectorI].name();
898 
899  if (patchSource.found("patchFields"))
900  {
901  dictionary patchFieldsDict = patchSource.subDict
902  (
903  "patchFields"
904  );
905 
906  if (sameGroup)
907  {
908  // Add coupleGroup to all entries
909  forAllIter(dictionary, patchFieldsDict, iter)
910  {
911  if (iter().isDict())
912  {
913  dictionary& dict = iter().dict();
914  dict.set("coupleGroup", groupName);
915  }
916  }
917 
918  const labelList& patchIDs =
919  pbm.groupPatchIDs()[groupName];
920 
921  forAll(patchIDs, i)
922  {
923  fvMeshTools::setPatchFields
924  (
925  mesh,
926  patchIDs[i],
927  patchFieldsDict
928  );
929  }
930  }
931  else
932  {
933  const word masterPatchName(groupName + "_master");
934  const word slavePatchName(groupName + "_slave");
935 
936  label patchiMaster = pbm.findPatchID(masterPatchName);
937  label patchiSlave = pbm.findPatchID(slavePatchName);
938 
939  fvMeshTools::setPatchFields
940  (
941  mesh,
942  patchiMaster,
943  patchFieldsDict
944  );
945 
946  fvMeshTools::setPatchFields
947  (
948  mesh,
949  patchiSlave,
950  patchFieldsDict
951  );
952  }
953  }
954  }
955  }
956  }
957 
958 
959  // Move mesh (since morphing might not do this)
960  if (map().hasMotionPoints())
961  {
962  mesh.movePoints(map().preMotionPoints());
963  }
964 
965 
966  // Remove any now zero-sized patches
967  filterPatches(mesh, bafflePatches);
968 
969 
970  if (overwrite)
971  {
972  mesh.setInstance(oldInstance);
973  }
974 
975  Info<< "Writing mesh to " << runTime.timeName() << endl;
976 
977  mesh.write();
978 
979  Info<< "End\n" << endl;
980 
981  return 0;
982 }
983 
984 
985 // ************************************************************************* //
Foam::surfaceFields.
A HashTable with keys but without contents.
Definition: HashSet.H:59
dictionary dict
static void reorderPatches(fvMesh &, const labelList &oldToNew, const label nPatches, const bool validBoundary)
Reorder and remove trailing patches. If validBoundary call is parallel.
Definition: fvMeshTools.C:324
wordList ReadFields(const Mesh &mesh, const IOobjectList &objects, PtrList< GeoField > &fields, const bool syncPar=true)
Read all fields of the specified type.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
void clearAddressing()
Clear addressing.
Definition: ZoneMesh.C:387
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
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:223
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
bool checkParallelSync(const bool report=false) const
Check whether all procs have all patches and in same order. Return.
const boolList & flipMap() const
Return face flip map.
Definition: faceZone.H:253
Class describing modification of a face.
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
label index() const
Return the index of this zone in zone list.
Definition: zone.H:161
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:453
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:633
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
Field reading functions for post-processing utilities.
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:116
Class holds all the necessary information for mapping fields associated with fvMesh.
Definition: fvMeshMapper.H:55
bool checkParallelSync(const bool report=false) const
Check whether all procs have all zones and in same order. Return.
Definition: ZoneMesh.C:428
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
A face addition data class. A face can be inflated either from a point or from another face and can e...
Definition: polyAddFace.H:51
bool inGroup(const word &) const
Test if in group.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
static void zeroPatchFields(fvMesh &mesh, const label patchi)
Change patchField to zero on registered fields.
Definition: fvMeshTools.C:238
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
Definition: faceZone.C:305
const word & name() const
Return name.
Definition: zone.H:150
label findZoneID(const word &zoneName) const
Find zone index given a name.
Definition: ZoneMesh.C:341
virtual label sizeBeforeMapping() const
Return size of field before mapping.
virtual void updateMesh(const mapPolyMesh &mpm)
Update mesh corresponding to the given map.
Definition: fvMesh.C:800
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
Definition: polyPatch.H:310
A class for handling words, derived from string.
Definition: word.H:59
static label addPatch(fvMesh &mesh, const polyPatch &patch, const dictionary &patchFieldDict, const word &defaultPatchFieldType, const bool validBoundary)
Add patch. Inserts patch before all processor patches.
Definition: fvMeshTools.C:32
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:97
bool set(const label) const
Is element set.
Definition: PtrListI.H:65
This boundary condition is not designed to be evaluated; it is assmued that the value is assigned via...
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:292
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
const word & name() const
Return name.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
bool found(const Key &) const
Return true if hashedEntry is found in table.
Foam::polyBoundaryMesh.
static const char nl
Definition: Ostream.H:262
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:237
const HashTable< labelList, word > & groupPatchIDs() const
Per patch group the patch indices.
IOobject dictIO(dictName, runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE)
void setInstance(const fileName &)
Set the instance for mesh files.
Definition: polyMeshIO.C:32
static autoPtr< faceSelection > New(const word &name, const fvMesh &mesh, const dictionary &dict)
Return a reference to the selected faceSelection.
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
Definition: fvMesh.C:721
label nFaces() const
A bit-packed bool list.
label patchi
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1023
#define WarningInFunction
Report a warning using Foam::Warning.
word dictName("noiseDict")
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:772
List< Key > sortedToc() const
Return the table of contents as a sorted list.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Direct mesh changes based on v1.3 polyTopoChange syntax.
void set(entry *)
Assign a new entry, overwrite any existing entry.
Definition: dictionary.C:864
messageStream Info
label n
const faceZoneMesh & faceZones() const
Return face zone mesh.
Definition: polyMesh.H:463
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:64
Mapping class for a fvPatchField.
Definition: fvPatchMapper.H:55
static void addNote(const string &)
Add extra notes for the usage information.
Definition: argList.C:124
autoPtr< mapPolyMesh > changeMesh(polyMesh &mesh, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
const wordList & inGroups() const
Return the optional groups patch belongs to.
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1017
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1004
Foam::argList args(argc, argv)
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:230
label findPatchID(const word &patchName) const
Find patch index given a name.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
label nInternalFaces() const
bool found
virtual bool write() const
Write mesh using IO settings from time.
Definition: fvMesh.C:870
const word & name() const
Return name.
Definition: IOobject.H:260
Namespace for OpenFOAM.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451