FreeFOAM The Cross-Platform CFD Toolkit
polyMesh.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 \*---------------------------------------------------------------------------*/
25 
26 #include "polyMesh.H"
27 #include <OpenFOAM/Time.H>
28 #include <OpenFOAM/cellIOList.H>
29 #include <OpenFOAM/SubList.H>
34 #include <OpenFOAM/OSspecific.H>
36 
37 #include <OpenFOAM/pointMesh.H>
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
42 
43 
46 
47 
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 
50 void Foam::polyMesh::calcDirections() const
51 {
52  for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
53  {
54  solutionD_[cmpt] = 1;
55  }
56 
57  // Knock out empty and wedge directions. Note:they will be present on all
58  // domains.
59 
60  label nEmptyPatches = 0;
61  label nWedgePatches = 0;
62 
63  vector emptyDirVec = vector::zero;
64  vector wedgeDirVec = vector::zero;
65 
67  {
68  if (boundaryMesh()[patchi].size())
69  {
70  if (isA<emptyPolyPatch>(boundaryMesh()[patchi]))
71  {
72  nEmptyPatches++;
73  emptyDirVec += sum(cmptMag(boundaryMesh()[patchi].faceAreas()));
74  }
75  else if (isA<wedgePolyPatch>(boundaryMesh()[patchi]))
76  {
77  const wedgePolyPatch& wpp = refCast<const wedgePolyPatch>
78  (
80  );
81 
82  nWedgePatches++;
83  wedgeDirVec += cmptMag(wpp.centreNormal());
84  }
85  }
86  }
87 
88  reduce(nEmptyPatches, maxOp<label>());
89  reduce(nWedgePatches, maxOp<label>());
90 
91  if (nEmptyPatches)
92  {
93  reduce(emptyDirVec, sumOp<vector>());
94 
95  emptyDirVec /= mag(emptyDirVec);
96 
97  for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
98  {
99  if (emptyDirVec[cmpt] > 1e-6)
100  {
101  solutionD_[cmpt] = -1;
102  }
103  else
104  {
105  solutionD_[cmpt] = 1;
106  }
107  }
108  }
109 
110 
111  // Knock out wedge directions
112 
113  geometricD_ = solutionD_;
114 
115  if (nWedgePatches)
116  {
117  reduce(wedgeDirVec, sumOp<vector>());
118 
119  wedgeDirVec /= mag(wedgeDirVec);
120 
121  for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
122  {
123  if (wedgeDirVec[cmpt] > 1e-6)
124  {
125  geometricD_[cmpt] = -1;
126  }
127  else
128  {
129  geometricD_[cmpt] = 1;
130  }
131  }
132  }
133 }
134 
135 
136 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
137 
138 Foam::polyMesh::polyMesh(const IOobject& io)
139 :
140  objectRegistry(io),
141  primitiveMesh(),
142  points_
143  (
144  IOobject
145  (
146  "points",
147  time().findInstance(meshDir(), "points"),
148  meshSubDir,
149  *this,
150  IOobject::MUST_READ,
151  IOobject::NO_WRITE
152  )
153  ),
154  faces_
155  (
156  IOobject
157  (
158  "faces",
159  time().findInstance(meshDir(), "faces"),
160  meshSubDir,
161  *this,
162  IOobject::MUST_READ,
163  IOobject::NO_WRITE
164  )
165  ),
166  owner_
167  (
168  IOobject
169  (
170  "owner",
171  time().findInstance(meshDir(), "faces"),
172  meshSubDir,
173  *this,
174  IOobject::READ_IF_PRESENT,
175  IOobject::NO_WRITE
176  )
177  ),
178  neighbour_
179  (
180  IOobject
181  (
182  "neighbour",
183  time().findInstance(meshDir(), "faces"),
184  meshSubDir,
185  *this,
186  IOobject::READ_IF_PRESENT,
187  IOobject::NO_WRITE
188  )
189  ),
190  clearedPrimitives_(false),
191  boundary_
192  (
193  IOobject
194  (
195  "boundary",
196  time().findInstance(meshDir(), "boundary"),
197  meshSubDir,
198  *this,
199  IOobject::MUST_READ,
200  IOobject::NO_WRITE
201  ),
202  *this
203  ),
204  bounds_(points_),
205  geometricD_(Vector<label>::zero),
206  solutionD_(Vector<label>::zero),
207  pointZones_
208  (
209  IOobject
210  (
211  "pointZones",
212  time().findInstance
213  (
214  meshDir(),
215  "pointZones",
216  IOobject::READ_IF_PRESENT
217  ),
218  meshSubDir,
219  *this,
220  IOobject::READ_IF_PRESENT,
221  IOobject::NO_WRITE
222  ),
223  *this
224  ),
225  faceZones_
226  (
227  IOobject
228  (
229  "faceZones",
230  time().findInstance
231  (
232  meshDir(),
233  "faceZones",
234  IOobject::READ_IF_PRESENT
235  ),
236  meshSubDir,
237  *this,
238  IOobject::READ_IF_PRESENT,
239  IOobject::NO_WRITE
240  ),
241  *this
242  ),
243  cellZones_
244  (
245  IOobject
246  (
247  "cellZones",
248  time().findInstance
249  (
250  meshDir(),
251  "cellZones",
252  IOobject::READ_IF_PRESENT
253  ),
254  meshSubDir,
255  *this,
256  IOobject::READ_IF_PRESENT,
257  IOobject::NO_WRITE
258  ),
259  *this
260  ),
261  globalMeshDataPtr_(NULL),
262  moving_(false),
263  changing_(false),
264  curMotionTimeIndex_(time().timeIndex()),
265  oldPointsPtr_(NULL)
266 {
267  if (exists(owner_.objectPath()))
268  {
269  initMesh();
270  }
271  else
272  {
273  cellIOList cLst
274  (
275  IOobject
276  (
277  "cells",
278  time().findInstance(meshDir(), "cells"),
279  meshSubDir,
280  *this,
283  )
284  );
285 
286  // Set the primitive mesh
287  initMesh(cLst);
288 
289  owner_.write();
290  neighbour_.write();
291  }
292 
293  // Calculate topology for the patches (processor-processor comms etc.)
294  boundary_.updateMesh();
295 
296  // Calculate the geometry for the patches (transformation tensors etc.)
297  boundary_.calcGeometry();
298 
299  // Warn if global empty mesh
300  if (returnReduce(nPoints(), sumOp<label>()) == 0)
301  {
302  WarningIn("polyMesh(const IOobject&)")
303  << "no points in mesh" << endl;
304  }
305  if (returnReduce(nCells(), sumOp<label>()) == 0)
306  {
307  WarningIn("polyMesh(const IOobject&)")
308  << "no cells in mesh" << endl;
309  }
310 }
311 
312 
313 Foam::polyMesh::polyMesh
314 (
315  const IOobject& io,
316  const Xfer<pointField>& points,
317  const Xfer<faceList>& faces,
318  const Xfer<labelList>& owner,
319  const Xfer<labelList>& neighbour,
320  const bool syncPar
321 )
322 :
323  objectRegistry(io),
324  primitiveMesh(),
325  points_
326  (
327  IOobject
328  (
329  "points",
330  instance(),
331  meshSubDir,
332  *this,
335  ),
336  points
337  ),
338  faces_
339  (
340  IOobject
341  (
342  "faces",
343  instance(),
344  meshSubDir,
345  *this,
348  ),
349  faces
350  ),
351  owner_
352  (
353  IOobject
354  (
355  "owner",
356  instance(),
357  meshSubDir,
358  *this,
361  ),
362  owner
363  ),
364  neighbour_
365  (
366  IOobject
367  (
368  "neighbour",
369  instance(),
370  meshSubDir,
371  *this,
374  ),
375  neighbour
376  ),
377  clearedPrimitives_(false),
378  boundary_
379  (
380  IOobject
381  (
382  "boundary",
383  instance(),
384  meshSubDir,
385  *this,
388  ),
389  *this,
390  0
391  ),
392  bounds_(points_, syncPar),
393  geometricD_(Vector<label>::zero),
394  solutionD_(Vector<label>::zero),
395  pointZones_
396  (
397  IOobject
398  (
399  "pointZones",
400  instance(),
401  meshSubDir,
402  *this,
405  ),
406  *this,
407  0
408  ),
409  faceZones_
410  (
411  IOobject
412  (
413  "faceZones",
414  instance(),
415  meshSubDir,
416  *this,
419  ),
420  *this,
421  0
422  ),
423  cellZones_
424  (
425  IOobject
426  (
427  "cellZones",
428  instance(),
429  meshSubDir,
430  *this,
433  ),
434  *this,
435  0
436  ),
437  globalMeshDataPtr_(NULL),
438  moving_(false),
439  changing_(false),
440  curMotionTimeIndex_(time().timeIndex()),
441  oldPointsPtr_(NULL)
442 {
443  // Check if the faces and cells are valid
444  forAll (faces_, faceI)
445  {
446  const face& curFace = faces_[faceI];
447 
448  if (min(curFace) < 0 || max(curFace) > points_.size())
449  {
451  (
452  "polyMesh::polyMesh\n"
453  "(\n"
454  " const IOobject& io,\n"
455  " const pointField& points,\n"
456  " const faceList& faces,\n"
457  " const cellList& cells\n"
458  ")\n"
459  ) << "Face " << faceI << "contains vertex labels out of range: "
460  << curFace << " Max point index = " << points_.size()
461  << abort(FatalError);
462  }
463  }
464 
465  // Set the primitive mesh
466  initMesh();
467 }
468 
469 
470 Foam::polyMesh::polyMesh
471 (
472  const IOobject& io,
473  const Xfer<pointField>& points,
474  const Xfer<faceList>& faces,
475  const Xfer<cellList>& cells,
476  const bool syncPar
477 )
478 :
479  objectRegistry(io),
480  primitiveMesh(),
481  points_
482  (
483  IOobject
484  (
485  "points",
486  instance(),
487  meshSubDir,
488  *this,
491  ),
492  points
493  ),
494  faces_
495  (
496  IOobject
497  (
498  "faces",
499  instance(),
500  meshSubDir,
501  *this,
504  ),
505  faces
506  ),
507  owner_
508  (
509  IOobject
510  (
511  "owner",
512  instance(),
513  meshSubDir,
514  *this,
517  ),
518  0
519  ),
520  neighbour_
521  (
522  IOobject
523  (
524  "neighbour",
525  instance(),
526  meshSubDir,
527  *this,
530  ),
531  0
532  ),
533  clearedPrimitives_(false),
534  boundary_
535  (
536  IOobject
537  (
538  "boundary",
539  instance(),
540  meshSubDir,
541  *this,
544  ),
545  *this,
546  0
547  ),
548  bounds_(points_, syncPar),
549  geometricD_(Vector<label>::zero),
550  solutionD_(Vector<label>::zero),
551  pointZones_
552  (
553  IOobject
554  (
555  "pointZones",
556  instance(),
557  meshSubDir,
558  *this,
561  ),
562  *this,
563  0
564  ),
565  faceZones_
566  (
567  IOobject
568  (
569  "faceZones",
570  instance(),
571  meshSubDir,
572  *this,
575  ),
576  *this,
577  0
578  ),
579  cellZones_
580  (
581  IOobject
582  (
583  "cellZones",
584  instance(),
585  meshSubDir,
586  *this,
589  ),
590  *this,
591  0
592  ),
593  globalMeshDataPtr_(NULL),
594  moving_(false),
595  changing_(false),
596  curMotionTimeIndex_(time().timeIndex()),
597  oldPointsPtr_(NULL)
598 {
599  // Check if faces are valid
600  forAll (faces_, faceI)
601  {
602  const face& curFace = faces_[faceI];
603 
604  if (min(curFace) < 0 || max(curFace) > points_.size())
605  {
607  (
608  "polyMesh::polyMesh\n"
609  "(\n"
610  " const IOobject&,\n"
611  " const Xfer<pointField>&,\n"
612  " const Xfer<faceList>&,\n"
613  " const Xfer<cellList>&\n"
614  ")\n"
615  ) << "Face " << faceI << "contains vertex labels out of range: "
616  << curFace << " Max point index = " << points_.size()
617  << abort(FatalError);
618  }
619  }
620 
621  // transfer in cell list
622  cellList cLst(cells);
623 
624  // Check if cells are valid
625  forAll (cLst, cellI)
626  {
627  const cell& curCell = cLst[cellI];
628 
629  if (min(curCell) < 0 || max(curCell) > faces_.size())
630  {
632  (
633  "polyMesh::polyMesh\n"
634  "(\n"
635  " const IOobject&,\n"
636  " const Xfer<pointField>&,\n"
637  " const Xfer<faceList>&,\n"
638  " const Xfer<cellList>&\n"
639  ")\n"
640  ) << "Cell " << cellI << "contains face labels out of range: "
641  << curCell << " Max face index = " << faces_.size()
642  << abort(FatalError);
643  }
644  }
645 
646  // Set the primitive mesh
647  initMesh(cLst);
648 }
649 
650 
652 (
653  const Xfer<pointField>& points,
654  const Xfer<faceList>& faces,
655  const Xfer<labelList>& owner,
656  const Xfer<labelList>& neighbour,
657  const labelList& patchSizes,
658  const labelList& patchStarts,
659  const bool validBoundary
660 )
661 {
662  // Clear addressing. Keep geometric props for mapping.
663  clearAddressing();
664 
665  // Take over new primitive data.
666  // Optimized to avoid overwriting data at all
667  if (&points)
668  {
669  points_.transfer(points());
670  bounds_ = boundBox(points_, validBoundary);
671  }
672 
673  if (&faces)
674  {
675  faces_.transfer(faces());
676  }
677 
678  if (&owner)
679  {
680  owner_.transfer(owner());
681  }
682 
683  if (&neighbour)
684  {
685  neighbour_.transfer(neighbour());
686  }
687 
688 
689  // Reset patch sizes and starts
690  forAll(boundary_, patchI)
691  {
692  boundary_[patchI] = polyPatch
693  (
694  boundary_[patchI].name(),
695  patchSizes[patchI],
696  patchStarts[patchI],
697  patchI,
698  boundary_
699  );
700  }
701 
702 
703  // Flags the mesh files as being changed
704  setInstance(time().timeName());
705 
706  // Check if the faces and cells are valid
707  forAll (faces_, faceI)
708  {
709  const face& curFace = faces_[faceI];
710 
711  if (min(curFace) < 0 || max(curFace) > points_.size())
712  {
714  (
715  "polyMesh::polyMesh::resetPrimitives\n"
716  "(\n"
717  " const Xfer<pointField>&,\n"
718  " const Xfer<faceList>&,\n"
719  " const Xfer<labelList>& owner,\n"
720  " const Xfer<labelList>& neighbour,\n"
721  " const labelList& patchSizes,\n"
722  " const labelList& patchStarts\n"
723  ")\n"
724  ) << "Face " << faceI << " contains vertex labels out of range: "
725  << curFace << " Max point index = " << points_.size()
726  << abort(FatalError);
727  }
728  }
729 
730 
731  // Set the primitive mesh from the owner_, neighbour_.
732  // Works out from patch end where the active faces stop.
733  initMesh();
734 
735 
736  if (validBoundary)
737  {
738  // Note that we assume that all the patches stay the same and are
739  // correct etc. so we can already use the patches to do
740  // processor-processor comms.
741 
742  // Calculate topology for the patches (processor-processor comms etc.)
743  boundary_.updateMesh();
744 
745  // Calculate the geometry for the patches (transformation tensors etc.)
746  boundary_.calcGeometry();
747 
748  // Warn if global empty mesh
749  if
750  (
751  (returnReduce(nPoints(), sumOp<label>()) == 0)
752  || (returnReduce(nCells(), sumOp<label>()) == 0)
753  )
754  {
756  (
757  "polyMesh::polyMesh::resetPrimitives\n"
758  "(\n"
759  " const Xfer<pointField>&,\n"
760  " const Xfer<faceList>&,\n"
761  " const Xfer<labelList>& owner,\n"
762  " const Xfer<labelList>& neighbour,\n"
763  " const labelList& patchSizes,\n"
764  " const labelList& patchStarts\n"
765  " const bool validBoundary\n"
766  ")\n"
767  ) << "no points or no cells in mesh" << endl;
768  }
769  }
770 }
771 
772 
773 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
774 
776 {
777  clearOut();
778  resetMotion();
779 }
780 
781 
782 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
783 
785 {
786  if (objectRegistry::dbDir() == defaultRegion)
787  {
788  return parent().dbDir();
789  }
790  else
791  {
792  return objectRegistry::dbDir();
793  }
794 }
795 
796 
798 {
799  return dbDir()/meshSubDir;
800 }
801 
802 
804 {
805  return points_.instance();
806 }
807 
808 
810 {
811  return faces_.instance();
812 }
813 
814 
816 {
817  if (geometricD_.x() == 0)
818  {
819  calcDirections();
820  }
821 
822  return geometricD_;
823 }
824 
825 
826 Foam::label Foam::polyMesh::nGeometricD() const
827 {
828  return cmptSum(geometricD() + Vector<label>::one)/2;
829 }
830 
831 
833 {
834  if (solutionD_.x() == 0)
835  {
836  calcDirections();
837  }
838 
839  return solutionD_;
840 }
841 
842 
843 Foam::label Foam::polyMesh::nSolutionD() const
844 {
845  return cmptSum(solutionD() + Vector<label>::one)/2;
846 }
847 
848 
849 // Add boundary patches. Constructor helper
851 (
852  const List<polyPatch*>& p,
853  const bool validBoundary
854 )
855 {
856  if (boundaryMesh().size())
857  {
859  (
860  "void polyMesh::addPatches(const List<polyPatch*>&, const bool)"
861  ) << "boundary already exists"
862  << abort(FatalError);
863  }
864 
865  // Reset valid directions
866  geometricD_ = Vector<label>::zero;
867  solutionD_ = Vector<label>::zero;
868 
869  boundary_.setSize(p.size());
870 
871  // Copy the patch pointers
872  forAll (p, pI)
873  {
874  boundary_.set(pI, p[pI]);
875  }
876 
877  // parallelData depends on the processorPatch ordering so force
878  // recalculation. Problem: should really be done in removeBoundary but
879  // there is some info in parallelData which might be interesting inbetween
880  // removeBoundary and addPatches.
881  deleteDemandDrivenData(globalMeshDataPtr_);
882 
883  if (validBoundary)
884  {
885  // Calculate topology for the patches (processor-processor comms etc.)
886  boundary_.updateMesh();
887 
888  // Calculate the geometry for the patches (transformation tensors etc.)
889  boundary_.calcGeometry();
890 
891  boundary_.checkDefinition();
892  }
893 }
894 
895 
896 // Add mesh zones. Constructor helper
898 (
899  const List<pointZone*>& pz,
900  const List<faceZone*>& fz,
901  const List<cellZone*>& cz
902 )
903 {
904  if (pointZones().size() || faceZones().size() || cellZones().size())
905  {
907  (
908  "void addZones\n"
909  "(\n"
910  " const List<pointZone*>&,\n"
911  " const List<faceZone*>&,\n"
912  " const List<cellZone*>&\n"
913  ")"
914  ) << "point, face or cell zone already exists"
915  << abort(FatalError);
916  }
917 
918  // Point zones
919  if (pz.size())
920  {
921  pointZones_.setSize(pz.size());
922 
923  // Copy the zone pointers
924  forAll (pz, pI)
925  {
926  pointZones_.set(pI, pz[pI]);
927  }
928 
929  pointZones_.writeOpt() = IOobject::AUTO_WRITE;
930  }
931 
932  // Face zones
933  if (fz.size())
934  {
935  faceZones_.setSize(fz.size());
936 
937  // Copy the zone pointers
938  forAll (fz, fI)
939  {
940  faceZones_.set(fI, fz[fI]);
941  }
942 
943  faceZones_.writeOpt() = IOobject::AUTO_WRITE;
944  }
945 
946  // Cell zones
947  if (cz.size())
948  {
949  cellZones_.setSize(cz.size());
950 
951  // Copy the zone pointers
952  forAll (cz, cI)
953  {
954  cellZones_.set(cI, cz[cI]);
955  }
956 
957  cellZones_.writeOpt() = IOobject::AUTO_WRITE;
958  }
959 }
960 
961 
963 {
964  if (clearedPrimitives_)
965  {
966  FatalErrorIn("const pointField& polyMesh::points() const")
967  << "points deallocated"
968  << abort(FatalError);
969  }
970 
971  return points_;
972 }
973 
974 
976 {
977  if (clearedPrimitives_)
978  {
979  FatalErrorIn("const faceList& polyMesh::faces() const")
980  << "faces deallocated"
981  << abort(FatalError);
982  }
983 
984  return faces_;
985 }
986 
987 
989 {
990  return owner_;
991 }
992 
993 
995 {
996  return neighbour_;
997 }
998 
999 
1000 // Return old mesh motion points
1002 {
1003  if (!oldPointsPtr_)
1004  {
1005  if (debug)
1006  {
1007  WarningIn("const pointField& polyMesh::oldPoints() const")
1008  << "Old points not available. Forcing storage of old points"
1009  << endl;
1010  }
1011 
1012  oldPointsPtr_ = new pointField(points_);
1013  curMotionTimeIndex_ = time().timeIndex();
1014  }
1015 
1016  return *oldPointsPtr_;
1017 }
1018 
1019 
1022  const pointField& newPoints
1023 )
1024 {
1025  if (debug)
1026  {
1027  Info<< "tmp<scalarField> polyMesh::movePoints(const pointField&) : "
1028  << " Moving points for time " << time().value()
1029  << " index " << time().timeIndex() << endl;
1030  }
1031 
1032  moving(true);
1033 
1034  // Pick up old points
1035  if (curMotionTimeIndex_ != time().timeIndex())
1036  {
1037  // Mesh motion in the new time step
1038  deleteDemandDrivenData(oldPointsPtr_);
1039  oldPointsPtr_ = new pointField(points_);
1040  curMotionTimeIndex_ = time().timeIndex();
1041  }
1042 
1043  points_ = newPoints;
1044 
1045  if (debug)
1046  {
1047  // Check mesh motion
1048  if (primitiveMesh::checkMeshMotion(points_, true))
1049  {
1050  Info<< "tmp<scalarField> polyMesh::movePoints"
1051  << "(const pointField&) : "
1052  << "Moving the mesh with given points will "
1053  << "invalidate the mesh." << nl
1054  << "Mesh motion should not be executed." << endl;
1055  }
1056  }
1057 
1058  points_.writeOpt() = IOobject::AUTO_WRITE;
1059  points_.instance() = time().timeName();
1060 
1061 
1063  (
1064  points_,
1065  oldPoints()
1066  );
1067 
1068  // Adjust parallel shared points
1069  if (globalMeshDataPtr_)
1070  {
1071  globalMeshDataPtr_->movePoints(points_);
1072  }
1073 
1074  // Force recalculation of all geometric data with new points
1075 
1076  bounds_ = boundBox(points_);
1077  boundary_.movePoints(points_);
1078 
1079  pointZones_.movePoints(points_);
1080  faceZones_.movePoints(points_);
1081  cellZones_.movePoints(points_);
1082 
1083  // Reset valid directions (could change with rotation)
1084  geometricD_ = Vector<label>::zero;
1085  solutionD_ = Vector<label>::zero;
1086 
1087 
1088  // Hack until proper callbacks. Below are all the polyMeh MeshObjects with a
1089  // movePoints function.
1090 
1091  // pointMesh
1092  if (thisDb().foundObject<pointMesh>(pointMesh::typeName))
1093  {
1094  const_cast<pointMesh&>
1095  (
1096  thisDb().lookupObject<pointMesh>
1097  (
1098  pointMesh::typeName
1099  )
1100  ).movePoints(points_);
1101  }
1102 
1103  return sweptVols;
1104 }
1105 
1106 
1107 // Reset motion by deleting old points
1109 {
1110  curMotionTimeIndex_ = 0;
1111  deleteDemandDrivenData(oldPointsPtr_);
1112 }
1113 
1114 
1115 // Return parallel info
1117 {
1118  if (!globalMeshDataPtr_)
1119  {
1120  if (debug)
1121  {
1122  Pout<< "polyMesh::globalData() const : "
1123  << "Constructing parallelData from processor topology" << nl
1124  << "This needs the patch faces to be correctly matched"
1125  << endl;
1126  }
1127  // Construct globalMeshData using processorPatch information only.
1128  globalMeshDataPtr_ = new globalMeshData(*this);
1129  }
1130 
1131  return *globalMeshDataPtr_;
1132 }
1133 
1134 
1135 // Remove all files and some subdirs (eg, sets)
1136 void Foam::polyMesh::removeFiles(const fileName& instanceDir) const
1137 {
1138  fileName meshFilesPath = thisDb().path()/instanceDir/meshDir();
1139 
1140  rm(meshFilesPath/"points");
1141  rm(meshFilesPath/"faces");
1142  rm(meshFilesPath/"owner");
1143  rm(meshFilesPath/"neighbour");
1144  rm(meshFilesPath/"cells");
1145  rm(meshFilesPath/"boundary");
1146  rm(meshFilesPath/"pointZones");
1147  rm(meshFilesPath/"faceZones");
1148  rm(meshFilesPath/"cellZones");
1149  rm(meshFilesPath/"meshModifiers");
1150  rm(meshFilesPath/"parallelData");
1151 
1152  // remove subdirectories
1153  if (isDir(meshFilesPath/"sets"))
1154  {
1155  rmDir(meshFilesPath/"sets");
1156  }
1157 }
1158 
1160 {
1161  removeFiles(instance());
1162 }
1163 
1164 
1165 // ************************ vim: set sw=4 sts=4 et: ************************ //