FreeFOAM The Cross-Platform CFD Toolkit
createPatch.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 Application
25  createPatch
26 
27 Description
28  Utility to create patches out of selected boundary faces.
29 
30  Faces come either from existing patches or from a faceSet. More
31  specifically it:
32  - creates new patches (from selected boundary faces). Synchronise faces
33  on coupled patches.
34  - synchronises points on coupled boundaries
35  - remove patches with 0 faces in them
36 
37 Usage
38 
39  - createPatch [OPTIONS]
40 
41  @param -overwrite \n
42  Overwrite existing data.
43 
44  @param -case <dir>\n
45  Case directory.
46 
47  @param -parallel \n
48  Run in parallel.
49 
50  @param -help \n
51  Display help message.
52 
53  @param -doc \n
54  Display Doxygen API documentation page for this application.
55 
56  @param -srcDoc \n
57  Display Doxygen source documentation page for this application.
58 
59 \*---------------------------------------------------------------------------*/
60 
62 #include <OpenFOAM/syncTools.H>
63 #include <OpenFOAM/argList.H>
64 #include <OpenFOAM/polyMesh.H>
65 #include <OpenFOAM/Time.H>
66 #include <OpenFOAM/SortableList.H>
67 #include <OpenFOAM/OFstream.H>
68 #include <meshTools/meshTools.H>
69 #include <meshTools/faceSet.H>
70 #include <OpenFOAM/IOPtrList.H>
73 
74 using namespace Foam;
75 
76 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
77 
78 namespace Foam
79 {
81 }
82 
83 // Combine operator to synchronise points. We choose point nearest to origin so
84 // we can use e.g. great,great,great as null value.
85 class nearestEqOp
86 {
87 
88 public:
89 
90  void operator()(vector& x, const vector& y) const
91  {
92  if (magSqr(y) < magSqr(x))
93  {
94  x = y;
95  }
96  }
97 };
98 
99 
100 void changePatchID
101 (
102  const polyMesh& mesh,
103  const label faceID,
104  const label patchID,
105  polyTopoChange& meshMod
106 )
107 {
108  const label zoneID = mesh.faceZones().whichZone(faceID);
109 
110  bool zoneFlip = false;
111 
112  if (zoneID >= 0)
113  {
114  const faceZone& fZone = mesh.faceZones()[zoneID];
115 
116  zoneFlip = fZone.flipMap()[fZone.whichFace(faceID)];
117  }
118 
119  meshMod.setAction
120  (
122  (
123  mesh.faces()[faceID], // face
124  faceID, // face ID
125  mesh.faceOwner()[faceID], // owner
126  -1, // neighbour
127  false, // flip flux
128  patchID, // patch ID
129  false, // remove from zone
130  zoneID, // zone ID
131  zoneFlip // zone flip
132  )
133  );
134 }
135 
136 
137 // Filter out the empty patches.
138 void filterPatches(polyMesh& mesh)
139 {
140  const polyBoundaryMesh& patches = mesh.boundaryMesh();
141 
142  // Patches to keep
143  DynamicList<polyPatch*> allPatches(patches.size());
144 
145  label nOldPatches = returnReduce(patches.size(), sumOp<label>());
146 
147  // Copy old patches.
148  forAll(patches, patchI)
149  {
150  const polyPatch& pp = patches[patchI];
151 
152  // Note: reduce possible since non-proc patches guaranteed in same order
153  if (!isA<processorPolyPatch>(pp))
154  {
155  if (returnReduce(pp.size(), sumOp<label>()) > 0)
156  {
157  allPatches.append
158  (
159  pp.clone
160  (
161  patches,
162  allPatches.size(),
163  pp.size(),
164  pp.start()
165  ).ptr()
166  );
167  }
168  else
169  {
170  Info<< "Removing empty patch " << pp.name() << " at position "
171  << patchI << endl;
172  }
173  }
174  }
175  // Copy non-empty processor patches
176  forAll(patches, patchI)
177  {
178  const polyPatch& pp = patches[patchI];
179 
180  if (isA<processorPolyPatch>(pp))
181  {
182  if (pp.size())
183  {
184  allPatches.append
185  (
186  pp.clone
187  (
188  patches,
189  allPatches.size(),
190  pp.size(),
191  pp.start()
192  ).ptr()
193  );
194  }
195  else
196  {
197  Info<< "Removing empty processor patch " << pp.name()
198  << " at position " << patchI << endl;
199  }
200  }
201  }
202 
203  label nAllPatches = returnReduce(allPatches.size(), sumOp<label>());
204  if (nAllPatches != nOldPatches)
205  {
206  Info<< "Removing patches." << endl;
207  allPatches.shrink();
208  mesh.removeBoundary();
209  mesh.addPatches(allPatches);
210  }
211  else
212  {
213  Info<< "No patches removed." << endl;
214  }
215 }
216 
217 
218 // Dump for all patches the current match
219 void dumpCyclicMatch(const fileName& prefix, const polyMesh& mesh)
220 {
221  const polyBoundaryMesh& patches = mesh.boundaryMesh();
222 
223  forAll(patches, patchI)
224  {
225  if (isA<cyclicPolyPatch>(patches[patchI]))
226  {
227  const cyclicPolyPatch& cycPatch =
228  refCast<const cyclicPolyPatch>(patches[patchI]);
229 
230  label halfSize = cycPatch.size()/2;
231 
232  // Dump halves
233  {
234  OFstream str(prefix+cycPatch.name()+"_half0.obj");
235  Pout<< "Dumping " << cycPatch.name()
236  << " half0 faces to " << str.name() << endl;
238  (
239  str,
240  static_cast<faceList>
241  (
243  (
244  cycPatch,
245  halfSize
246  )
247  ),
248  cycPatch.points()
249  );
250  }
251  {
252  OFstream str(prefix+cycPatch.name()+"_half1.obj");
253  Pout<< "Dumping " << cycPatch.name()
254  << " half1 faces to " << str.name() << endl;
256  (
257  str,
258  static_cast<faceList>
259  (
261  (
262  cycPatch,
263  halfSize,
264  halfSize
265  )
266  ),
267  cycPatch.points()
268  );
269  }
270 
271 
272  // Lines between corresponding face centres
273  OFstream str(prefix+cycPatch.name()+"_match.obj");
274  label vertI = 0;
275 
276  Pout<< "Dumping cyclic match as lines between face centres to "
277  << str.name() << endl;
278 
279  for (label faceI = 0; faceI < halfSize; faceI++)
280  {
281  const point& fc0 = mesh.faceCentres()[cycPatch.start()+faceI];
282  meshTools::writeOBJ(str, fc0);
283  vertI++;
284 
285  label nbrFaceI = halfSize + faceI;
286  const point& fc1 =
287  mesh.faceCentres()[cycPatch.start()+nbrFaceI];
288  meshTools::writeOBJ(str, fc1);
289  vertI++;
290 
291  str<< "l " << vertI-1 << ' ' << vertI << nl;
292  }
293  }
294  }
295 }
296 
297 
298 void separateList
299 (
300  const vectorField& separation,
301  UList<vector>& field
302 )
303 {
304  if (separation.size() == 1)
305  {
306  // Single value for all.
307 
308  forAll(field, i)
309  {
310  field[i] += separation[0];
311  }
312  }
313  else if (separation.size() == field.size())
314  {
315  forAll(field, i)
316  {
317  field[i] += separation[i];
318  }
319  }
320  else
321  {
323  (
324  "separateList(const vectorField&, UList<vector>&)"
325  ) << "Sizes of field and transformation not equal. field:"
326  << field.size() << " transformation:" << separation.size()
327  << abort(FatalError);
328  }
329 }
330 
331 
332 // Synchronise points on both sides of coupled boundaries.
333 template <class CombineOp>
334 void syncPoints
335 (
336  const polyMesh& mesh,
338  const CombineOp& cop,
339  const point& nullValue
340 )
341 {
342  if (points.size() != mesh.nPoints())
343  {
345  (
346  "syncPoints"
347  "(const polyMesh&, pointField&, const CombineOp&, const point&)"
348  ) << "Number of values " << points.size()
349  << " is not equal to the number of points in the mesh "
350  << mesh.nPoints() << abort(FatalError);
351  }
352 
353  const polyBoundaryMesh& patches = mesh.boundaryMesh();
354 
355  // Is there any coupled patch with transformation?
356  bool hasTransformation = false;
357 
358  if (Pstream::parRun())
359  {
360  // Send
361 
362  forAll(patches, patchI)
363  {
364  const polyPatch& pp = patches[patchI];
365 
366  if
367  (
368  isA<processorPolyPatch>(pp)
369  && pp.nPoints() > 0
370  && refCast<const processorPolyPatch>(pp).owner()
371  )
372  {
373  const processorPolyPatch& procPatch =
374  refCast<const processorPolyPatch>(pp);
375 
376  // Get data per patchPoint in neighbouring point numbers.
377  pointField patchInfo(procPatch.nPoints(), nullValue);
378 
379  const labelList& meshPts = procPatch.meshPoints();
380  const labelList& nbrPts = procPatch.neighbPoints();
381 
382  forAll(nbrPts, pointI)
383  {
384  label nbrPointI = nbrPts[pointI];
385  if (nbrPointI >= 0 && nbrPointI < patchInfo.size())
386  {
387  patchInfo[nbrPointI] = points[meshPts[pointI]];
388  }
389  }
390 
391  OPstream toNbr(Pstream::blocking, procPatch.neighbProcNo());
392  toNbr << patchInfo;
393  }
394  }
395 
396 
397  // Receive and set.
398 
399  forAll(patches, patchI)
400  {
401  const polyPatch& pp = patches[patchI];
402 
403  if
404  (
405  isA<processorPolyPatch>(pp)
406  && pp.nPoints() > 0
407  && !refCast<const processorPolyPatch>(pp).owner()
408  )
409  {
410  const processorPolyPatch& procPatch =
411  refCast<const processorPolyPatch>(pp);
412 
413  pointField nbrPatchInfo(procPatch.nPoints());
414  {
415  // We do not know the number of points on the other side
416  // so cannot use Pstream::read.
417  IPstream fromNbr
418  (
420  procPatch.neighbProcNo()
421  );
422  fromNbr >> nbrPatchInfo;
423  }
424  // Null any value which is not on neighbouring processor
425  nbrPatchInfo.setSize(procPatch.nPoints(), nullValue);
426 
427  if (!procPatch.parallel())
428  {
429  hasTransformation = true;
430  transformList(procPatch.forwardT(), nbrPatchInfo);
431  }
432  else if (procPatch.separated())
433  {
434  hasTransformation = true;
435  separateList(-procPatch.separation(), nbrPatchInfo);
436  }
437 
438  const labelList& meshPts = procPatch.meshPoints();
439 
440  forAll(meshPts, pointI)
441  {
442  label meshPointI = meshPts[pointI];
443  points[meshPointI] = nbrPatchInfo[pointI];
444  }
445  }
446  }
447  }
448 
449  // Do the cyclics.
450  forAll(patches, patchI)
451  {
452  const polyPatch& pp = patches[patchI];
453 
454  if (isA<cyclicPolyPatch>(pp))
455  {
456  const cyclicPolyPatch& cycPatch =
457  refCast<const cyclicPolyPatch>(pp);
458 
459  const edgeList& coupledPoints = cycPatch.coupledPoints();
460  const labelList& meshPts = cycPatch.meshPoints();
461 
462  pointField half0Values(coupledPoints.size());
463 
464  forAll(coupledPoints, i)
465  {
466  const edge& e = coupledPoints[i];
467  label point0 = meshPts[e[0]];
468  half0Values[i] = points[point0];
469  }
470 
471  if (!cycPatch.parallel())
472  {
473  hasTransformation = true;
474  transformList(cycPatch.reverseT(), half0Values);
475  }
476  else if (cycPatch.separated())
477  {
478  hasTransformation = true;
479  const vectorField& v = cycPatch.coupledPolyPatch::separation();
480  separateList(v, half0Values);
481  }
482 
483  forAll(coupledPoints, i)
484  {
485  const edge& e = coupledPoints[i];
486  label point1 = meshPts[e[1]];
487  points[point1] = half0Values[i];
488  }
489  }
490  }
491 
492  //- Note: hasTransformation is only used for warning messages so
493  // reduction not strictly nessecary.
494  //reduce(hasTransformation, orOp<bool>());
495 
496  // Synchronize multiple shared points.
497  const globalMeshData& pd = mesh.globalData();
498 
499  if (pd.nGlobalPoints() > 0)
500  {
501  if (hasTransformation)
502  {
503  WarningIn
504  (
505  "syncPoints"
506  "(const polyMesh&, pointField&, const CombineOp&, const point&)"
507  ) << "There are decomposed cyclics in this mesh with"
508  << " transformations." << endl
509  << "This is not supported. The result will be incorrect"
510  << endl;
511  }
512 
513 
514  // Values on shared points.
515  pointField sharedPts(pd.nGlobalPoints(), nullValue);
516 
517  forAll(pd.sharedPointLabels(), i)
518  {
519  label meshPointI = pd.sharedPointLabels()[i];
520  // Fill my entries in the shared points
521  sharedPts[pd.sharedPointAddr()[i]] = points[meshPointI];
522  }
523 
524  // Combine on master.
525  Pstream::listCombineGather(sharedPts, cop);
526  Pstream::listCombineScatter(sharedPts);
527 
528  // Now we will all have the same information. Merge it back with
529  // my local information.
530  forAll(pd.sharedPointLabels(), i)
531  {
532  label meshPointI = pd.sharedPointLabels()[i];
533  points[meshPointI] = sharedPts[pd.sharedPointAddr()[i]];
534  }
535  }
536 }
537 
538 
539 // Main program:
540 
541 int main(int argc, char *argv[])
542 {
543 # include <OpenFOAM/addRegionOption.H>
544  argList::validOptions.insert("overwrite", "");
545 
546 # include <OpenFOAM/setRootCase.H>
547 # include <OpenFOAM/createTime.H>
548  runTime.functionObjects().off();
549 
550  Foam::word meshRegionName = polyMesh::defaultRegion;
551  args.optionReadIfPresent("region", meshRegionName);
552 
553  const bool overwrite = args.optionFound("overwrite");
554 
555  Info<< "Reading createPatchDict." << nl << endl;
556 
557  IOdictionary dict
558  (
559  IOobject
560  (
561  "createPatchDict",
562  runTime.system(),
563  (
564  meshRegionName != polyMesh::defaultRegion
565  ? meshRegionName
566  : word::null
567  ),
568  runTime,
571  false
572  )
573  );
574 
575 
576  // Whether to synchronise points
577  const Switch pointSync(dict.lookup("pointSync"));
578 
579 
580  // Set the matching tolerance so we can read illegal meshes
581  scalar tol = readScalar(dict.lookup("matchTolerance"));
582  Info<< "Using relative tolerance " << tol
583  << " to match up faces and points" << nl << endl;
585 
587 
588  const word oldInstance = mesh.pointsInstance();
589 
590  const polyBoundaryMesh& patches = mesh.boundaryMesh();
591 
592  // If running parallel check same patches everywhere
593  patches.checkParallelSync(true);
594 
595 
596  dumpCyclicMatch("initial_", mesh);
597 
598  // Read patch construct info from dictionary
599  PtrList<dictionary> patchSources(dict.lookup("patchInfo"));
600 
601 
602 
603  // 1. Add all new patches
604  // ~~~~~~~~~~~~~~~~~~~~~~
605 
606  if (patchSources.size())
607  {
608  // Old and new patches.
609  DynamicList<polyPatch*> allPatches(patches.size()+patchSources.size());
610 
611  label startFaceI = mesh.nInternalFaces();
612 
613  // Copy old patches.
614  forAll(patches, patchI)
615  {
616  const polyPatch& pp = patches[patchI];
617 
618  if (!isA<processorPolyPatch>(pp))
619  {
620  allPatches.append
621  (
622  pp.clone
623  (
624  patches,
625  patchI,
626  pp.size(),
627  startFaceI
628  ).ptr()
629  );
630  startFaceI += pp.size();
631  }
632  }
633 
634  forAll(patchSources, addedI)
635  {
636  const dictionary& dict = patchSources[addedI];
637 
638  word patchName(dict.lookup("name"));
639 
640  label destPatchI = patches.findPatchID(patchName);
641 
642  if (destPatchI == -1)
643  {
644  dictionary patchDict(dict.subDict("dictionary"));
645 
646  destPatchI = allPatches.size();
647 
648  Info<< "Adding new patch " << patchName
649  << " as patch " << destPatchI
650  << " from " << patchDict << endl;
651 
652  patchDict.set("nFaces", 0);
653  patchDict.set("startFace", startFaceI);
654 
655  // Add an empty patch.
656  allPatches.append
657  (
659  (
660  patchName,
661  patchDict,
662  destPatchI,
663  patches
664  ).ptr()
665  );
666  }
667  }
668 
669  // Copy old patches.
670  forAll(patches, patchI)
671  {
672  const polyPatch& pp = patches[patchI];
673 
674  if (isA<processorPolyPatch>(pp))
675  {
676  allPatches.append
677  (
678  pp.clone
679  (
680  patches,
681  patchI,
682  pp.size(),
683  startFaceI
684  ).ptr()
685  );
686  startFaceI += pp.size();
687  }
688  }
689 
690  allPatches.shrink();
691  mesh.removeBoundary();
692  mesh.addPatches(allPatches);
693 
694  Info<< endl;
695  }
696 
697 
698 
699  // 2. Repatch faces
700  // ~~~~~~~~~~~~~~~~
701 
702  polyTopoChange meshMod(mesh);
703 
704 
705  forAll(patchSources, addedI)
706  {
707  const dictionary& dict = patchSources[addedI];
708 
709  word patchName(dict.lookup("name"));
710 
711  label destPatchI = patches.findPatchID(patchName);
712 
713  if (destPatchI == -1)
714  {
715  FatalErrorIn(args.executable()) << "patch " << patchName
716  << " not added. Problem." << abort(FatalError);
717  }
718 
719  word sourceType(dict.lookup("constructFrom"));
720 
721  if (sourceType == "patches")
722  {
723  labelHashSet patchSources(patches.patchSet(dict.lookup("patches")));
724 
725  // Repatch faces of the patches.
726  forAllConstIter(labelHashSet, patchSources, iter)
727  {
728  const polyPatch& pp = patches[iter.key()];
729 
730  Info<< "Moving faces from patch " << pp.name()
731  << " to patch " << destPatchI << endl;
732 
733  forAll(pp, i)
734  {
735  changePatchID
736  (
737  mesh,
738  pp.start() + i,
739  destPatchI,
740  meshMod
741  );
742  }
743  }
744  }
745  else if (sourceType == "set")
746  {
747  word setName(dict.lookup("set"));
748 
749  faceSet faces(mesh, setName);
750 
751  Info<< "Read " << returnReduce(faces.size(), sumOp<label>())
752  << " faces from faceSet " << faces.name() << endl;
753 
754  // Sort (since faceSet contains faces in arbitrary order)
755  labelList faceLabels(faces.toc());
756 
757  SortableList<label> patchFaces(faceLabels);
758 
759  forAll(patchFaces, i)
760  {
761  label faceI = patchFaces[i];
762 
763  if (mesh.isInternalFace(faceI))
764  {
766  << "Face " << faceI << " specified in set "
767  << faces.name()
768  << " is not an external face of the mesh." << endl
769  << "This application can only repatch existing boundary"
770  << " faces." << exit(FatalError);
771  }
772 
773  changePatchID
774  (
775  mesh,
776  faceI,
777  destPatchI,
778  meshMod
779  );
780  }
781  }
782  else
783  {
785  << "Invalid source type " << sourceType << endl
786  << "Valid source types are 'patches' 'set'" << exit(FatalError);
787  }
788  }
789  Info<< endl;
790 
791 
792  // Change mesh, use inflation to reforce calculation of transformation
793  // tensors.
794  Info<< "Doing topology modification to order faces." << nl << endl;
795  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, true);
796  mesh.movePoints(map().preMotionPoints());
797 
798  dumpCyclicMatch("coupled_", mesh);
799 
800  // Synchronise points.
801  if (!pointSync)
802  {
803  Info<< "Not synchronising points." << nl << endl;
804  }
805  else
806  {
807  Info<< "Synchronising points." << nl << endl;
808 
809  // This is a bit tricky. Both normal and position might be out and
810  // current separation also includes the normal
811  // ( separation_ = (nf&(Cr - Cf))*nf ).
812 
813  // For processor patches:
814  // - disallow multiple separation/transformation. This basically
815  // excludes decomposed cyclics. Use the (probably 0) separation
816  // to align the points.
817  // For cyclic patches:
818  // - for separated ones use our own recalculated offset vector
819  // - for rotational ones use current one.
820 
821  forAll(mesh.boundaryMesh(), patchI)
822  {
823  const polyPatch& pp = mesh.boundaryMesh()[patchI];
824 
825  if (pp.size() && isA<coupledPolyPatch>(pp))
826  {
827  const coupledPolyPatch& cpp =
828  refCast<const coupledPolyPatch>(pp);
829 
830  if (cpp.separated())
831  {
832  Info<< "On coupled patch " << pp.name()
833  << " separation[0] was "
834  << cpp.separation()[0] << endl;
835 
836  if (isA<cyclicPolyPatch>(pp))
837  {
838  const cyclicPolyPatch& cycpp =
839  refCast<const cyclicPolyPatch>(pp);
840 
842  {
843  Info<< "On cyclic translation patch " << pp.name()
844  << " forcing uniform separation of "
845  << cycpp.separationVector() << endl;
846  const_cast<vectorField&>(cpp.separation()) =
847  pointField(1, cycpp.separationVector());
848  }
849  else
850  {
851  const_cast<vectorField&>(cpp.separation()) =
852  pointField
853  (
854  1,
855  pp[pp.size()/2].centre(mesh.points())
856  - pp[0].centre(mesh.points())
857  );
858  }
859  }
860  else
861  {
862  const_cast<vectorField&>(cpp.separation())
863  .setSize(1);
864  }
865  Info<< "On coupled patch " << pp.name()
866  << " forcing uniform separation of "
867  << cpp.separation() << endl;
868  }
869  else if (!cpp.parallel())
870  {
871  Info<< "On coupled patch " << pp.name()
872  << " forcing uniform rotation of "
873  << cpp.forwardT()[0] << endl;
874 
875  const_cast<tensorField&>
876  (
877  cpp.forwardT()
878  ).setSize(1);
879  const_cast<tensorField&>
880  (
881  cpp.reverseT()
882  ).setSize(1);
883 
884  Info<< "On coupled patch " << pp.name()
885  << " forcing uniform rotation of "
886  << cpp.forwardT() << endl;
887  }
888  }
889  }
890 
891  Info<< "Synchronising points." << endl;
892 
893  pointField newPoints(mesh.points());
894 
895  syncPoints
896  (
897  mesh,
898  newPoints,
899  nearestEqOp(),
900  point(GREAT, GREAT, GREAT)
901  );
902 
903  scalarField diff(mag(newPoints-mesh.points()));
904  Info<< "Points changed by average:" << gAverage(diff)
905  << " max:" << gMax(diff) << nl << endl;
906 
907  mesh.movePoints(newPoints);
908  }
909 
910  // 3. Remove zeros-sized patches
911  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
912 
913  Info<< "Removing patches with no faces in them." << nl<< endl;
914  filterPatches(mesh);
915 
916 
917  dumpCyclicMatch("final_", mesh);
918 
919 
920  // Set the precision of the points data to 10
922 
923  if (!overwrite)
924  {
925  runTime++;
926  }
927  else
928  {
929  mesh.setInstance(oldInstance);
930  }
931 
932  // Write resulting mesh
933  Info<< "Writing repatched mesh to " << runTime.timeName() << nl << endl;
934  mesh.write();
935 
936  Info<< "End\n" << endl;
937 
938  return 0;
939 }
940 
941 
942 // ************************ vim: set sw=4 sts=4 et: ************************ //