fvMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2023 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 \*---------------------------------------------------------------------------*/
25 
26 #include "fvMesh.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 #include "pointFields.H"
30 #include "slicedVolFields.H"
31 #include "slicedSurfaceFields.H"
32 #include "SubField.H"
33 #include "demandDrivenData.H"
34 #include "fvMeshLduAddressing.H"
35 #include "fvMeshTopoChanger.H"
36 #include "fvMeshDistributor.H"
37 #include "fvMeshMover.H"
38 #include "fvMeshStitcher.H"
39 #include "nonConformalFvPatch.H"
42 #include "polyTopoChangeMap.H"
43 #include "MapFvFields.H"
44 #include "fvMeshMapper.H"
45 #include "pointMesh.H"
46 #include "pointMeshMapper.H"
47 #include "MapPointField.H"
48 #include "meshObjects.H"
49 #include "HashPtrTable.H"
50 #include "CompactListList.H"
51 
52 #include "fvcSurfaceIntegrate.H"
53 #include "fvcReconstruct.H"
54 #include "surfaceInterpolate.H"
55 
56 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
57 
58 namespace Foam
59 {
61 }
62 
63 
64 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
65 
66 void Foam::fvMesh::clearGeomNotOldVol()
67 {
68  if (debug)
69  {
70  Pout<< FUNCTION_NAME << "clearGeomNotOldVol" << endl;
71  }
72 
74  <
75  fvMesh,
76  GeometricMeshObject,
77  MoveableMeshObject
78  >(*this);
79 
81  <
82  lduMesh,
83  GeometricMeshObject,
84  MoveableMeshObject
85  >(*this);
86 
88  deleteDemandDrivenData(SfSlicePtr_);
89  deleteDemandDrivenData(SfPtr_);
90  deleteDemandDrivenData(magSfSlicePtr_);
91  deleteDemandDrivenData(magSfPtr_);
92  deleteDemandDrivenData(CSlicePtr_);
94  deleteDemandDrivenData(CfSlicePtr_);
95  deleteDemandDrivenData(CfPtr_);
96 }
97 
98 
99 void Foam::fvMesh::updateGeomNotOldVol()
100 {
101  bool haveV = (VPtr_ != nullptr);
102  bool haveSf = (SfSlicePtr_ != nullptr || SfPtr_ != nullptr);
103  bool haveMagSf = (magSfSlicePtr_ != nullptr || magSfPtr_ != nullptr);
104  bool haveCP = (CSlicePtr_ != nullptr || CPtr_ != nullptr);
105  bool haveCf = (CfSlicePtr_ != nullptr || CfPtr_ != nullptr);
106 
107  clearGeomNotOldVol();
108 
109  // Now recreate the fields
110  if (haveV)
111  {
112  (void)V();
113  }
114  if (haveSf)
115  {
116  (void)Sf();
117  }
118  if (haveMagSf)
119  {
120  (void)magSf();
121  }
122  if (haveCP)
123  {
124  (void)C();
125  }
126  if (haveCf)
127  {
128  (void)Cf();
129  }
130 }
131 
132 
133 void Foam::fvMesh::clearGeom()
134 {
135  if (debug)
136  {
137  Pout<< FUNCTION_NAME << "Clearing geometric data" << endl;
138  }
139 
140  clearGeomNotOldVol();
141 
142  deleteDemandDrivenData(phiPtr_);
143  deleteDemandDrivenData(V0Ptr_);
144  deleteDemandDrivenData(V00Ptr_);
145 }
146 
147 
148 void Foam::fvMesh::clearAddressing(const bool isMeshUpdate)
149 {
150  if (debug)
151  {
152  Pout<< FUNCTION_NAME << "isMeshUpdate: " << isMeshUpdate << endl;
153  }
154 
155  if (isMeshUpdate)
156  {
157  // Part of a mesh update. Keep meshObjects that have an topoChange
158  // callback
160  <
161  fvMesh,
162  TopologicalMeshObject,
163  UpdateableMeshObject
164  >
165  (
166  *this
167  );
169  <
170  lduMesh,
171  TopologicalMeshObject,
172  UpdateableMeshObject
173  >
174  (
175  *this
176  );
177  }
178  else
179  {
180  meshObjects::clear<fvMesh, TopologicalMeshObject>(*this);
181  meshObjects::clear<lduMesh, TopologicalMeshObject>(*this);
182  }
183 
184  deleteDemandDrivenData(lduPtr_);
185  deleteDemandDrivenData(polyFacesBfPtr_);
186  deleteDemandDrivenData(polyBFaceOffsetsPtr_);
187  deleteDemandDrivenData(polyBFaceOffsetPatchesPtr_);
188  deleteDemandDrivenData(polyBFaceOffsetPatchFacesPtr_);
189  deleteDemandDrivenData(polyBFacePatchesPtr_);
190  deleteDemandDrivenData(polyBFacePatchFacesPtr_);
191 }
192 
193 
194 void Foam::fvMesh::storeOldVol(const scalarField& V)
195 {
196  if (curTimeIndex_ < time().timeIndex())
197  {
198  if (debug)
199  {
201  << " Storing old time volumes since from time " << curTimeIndex_
202  << " and time now " << time().timeIndex()
203  << " V:" << V.size()
204  << endl;
205  }
206 
207  if (V00Ptr_ && V0Ptr_)
208  {
209  // Copy V0 into V00 storage
210  *V00Ptr_ = *V0Ptr_;
211  }
212 
213  if (V0Ptr_)
214  {
215  // Copy V into V0 storage
216  V0Ptr_->scalarField::operator=(V);
217  }
218  else
219  {
220  // Allocate V0 storage, fill with V
221  V0Ptr_ = new DimensionedField<scalar, volMesh>
222  (
223  IOobject
224  (
225  "V0",
226  time().name(),
227  *this,
230  false
231  ),
232  *this,
233  dimVolume
234  );
235  scalarField& V0 = *V0Ptr_;
236  // Note: V0 now sized with current mesh, not with (potentially
237  // different size) V.
238  V0.setSize(V.size());
239  V0 = V;
240  }
241 
242  if (debug)
243  {
245  << " Stored old time volumes V0:" << V0Ptr_->size()
246  << endl;
247  if (V00Ptr_)
248  {
250  << " Stored oldold time volumes V00:" << V00Ptr_->size()
251  << endl;
252  }
253  }
254  }
255 }
256 
257 
258 void Foam::fvMesh::storeOldTimeFields()
259 {
260  storeOldTimeFields<PointField>();
261  storeOldTimeFields<VolField>();
262  storeOldTimeFields<SurfaceField>();
263 }
264 
265 
266 void Foam::fvMesh::nullOldestTimeFields()
267 {
268  nullOldestTimeFields<PointField>();
269  nullOldestTimeFields<VolField>();
270  nullOldestTimeFields<SurfaceField>();
271 }
272 
273 
275 {
276  clearGeom();
277 
279 
280  clearAddressing();
281 
283 }
284 
285 
286 Foam::wordList Foam::fvMesh::polyFacesPatchTypes() const
287 {
288  wordList wantedPatchTypes
289  (
290  boundary().size(),
292  );
293 
295  {
296  const fvPatch& fvp = boundary()[patchi];
297 
298  if (isA<nonConformalFvPatch>(fvp))
299  {
300  wantedPatchTypes[patchi] =
302  }
303  }
304 
305  return wantedPatchTypes;
306 }
307 
308 
309 Foam::surfaceLabelField::Boundary& Foam::fvMesh::polyFacesBfRef()
310 {
311  if (!polyFacesBfPtr_)
312  {
313  polyFacesBf();
314  }
315 
316  return *polyFacesBfPtr_;
317 }
318 
319 
320 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
321 
323 (
324  const IOobject& io,
325  const bool changers,
326  const stitchType stitch
327 )
328 :
329  polyMesh(io),
330  surfaceInterpolation(*this),
331  data(static_cast<const objectRegistry&>(*this)),
332  boundary_(*this, boundaryMesh()),
333  stitcher_(fvMeshStitcher::New(*this, changers).ptr()),
334  topoChanger_(nullptr),
335  distributor_(nullptr),
336  mover_(nullptr),
337  lduPtr_(nullptr),
338  polyFacesBfPtr_(nullptr),
339  polyBFaceOffsetsPtr_(nullptr),
340  polyBFaceOffsetPatchesPtr_(nullptr),
341  polyBFaceOffsetPatchFacesPtr_(nullptr),
342  polyBFacePatchesPtr_(nullptr),
343  polyBFacePatchFacesPtr_(nullptr),
344  curTimeIndex_(time().timeIndex()),
345  VPtr_(nullptr),
346  V0Ptr_(nullptr),
347  V00Ptr_(nullptr),
348  SfSlicePtr_(nullptr),
349  SfPtr_(nullptr),
350  magSfSlicePtr_(nullptr),
351  magSfPtr_(nullptr),
352  CSlicePtr_(nullptr),
353  CPtr_(nullptr),
354  CfSlicePtr_(nullptr),
355  CfPtr_(nullptr),
356  phiPtr_(nullptr)
357 {
358  if (debug)
359  {
360  Pout<< FUNCTION_NAME << "Constructing fvMesh from IOobject" << endl;
361  }
362 
363  // Stitch or Re-stitch if necessary
364  if (stitch != stitchType::none)
365  {
366  stitcher_->connect(false, stitch == stitchType::geometric, true);
367  }
368 
369  // Construct changers
370  if (changers)
371  {
372  topoChanger_.set(fvMeshTopoChanger::New(*this).ptr());
373  distributor_.set(fvMeshDistributor::New(*this).ptr());
374  mover_.set(fvMeshMover::New(*this).ptr());
375 
376  // Check the existence of the cell volumes and read if present
377  // and set the storage of V00
378  if (fileHandler().isFile(time().timePath()/"V0"))
379  {
381  (
382  IOobject
383  (
384  "V0",
385  time().name(),
386  *this,
389  false
390  ),
391  *this
392  );
393 
394  V00();
395  }
396 
397  // Check the existence of the mesh fluxes and read if present
398  if (fileHandler().isFile(time().timePath()/"meshPhi"))
399  {
400  phiPtr_ = new surfaceScalarField
401  (
402  IOobject
403  (
404  "meshPhi",
405  time().name(),
406  *this,
409  true
410  ),
411  *this
412  );
413  }
414  }
415 }
416 
417 
419 (
420  const IOobject& io,
421  pointField&& points,
422  const cellShapeList& shapes,
423  const faceListList& boundaryFaces,
424  const wordList& boundaryPatchNames,
425  const PtrList<dictionary>& boundaryDicts,
426  const word& defaultBoundaryPatchName,
427  const word& defaultBoundaryPatchType,
428  const bool syncPar
429 )
430 :
431  polyMesh
432  (
433  io,
434  std::move(points),
435  shapes,
436  boundaryFaces,
437  boundaryPatchNames,
438  boundaryDicts,
439  defaultBoundaryPatchName,
440  defaultBoundaryPatchType,
441  syncPar
442  ),
443  surfaceInterpolation(*this),
444  data(static_cast<const objectRegistry&>(*this)),
445  boundary_(*this, boundaryMesh()),
446  stitcher_(nullptr),
447  topoChanger_(nullptr),
448  distributor_(nullptr),
449  mover_(nullptr),
450  lduPtr_(nullptr),
451  polyFacesBfPtr_(nullptr),
452  polyBFaceOffsetsPtr_(nullptr),
453  polyBFaceOffsetPatchesPtr_(nullptr),
454  polyBFaceOffsetPatchFacesPtr_(nullptr),
455  polyBFacePatchesPtr_(nullptr),
456  polyBFacePatchFacesPtr_(nullptr),
457  curTimeIndex_(time().timeIndex()),
458  VPtr_(nullptr),
459  V0Ptr_(nullptr),
460  V00Ptr_(nullptr),
461  SfSlicePtr_(nullptr),
462  SfPtr_(nullptr),
463  magSfSlicePtr_(nullptr),
464  magSfPtr_(nullptr),
465  CSlicePtr_(nullptr),
466  CPtr_(nullptr),
467  CfSlicePtr_(nullptr),
468  CfPtr_(nullptr),
469  phiPtr_(nullptr)
470 {
471  if (debug)
472  {
473  Pout<< FUNCTION_NAME << "Constructing fvMesh from cellShapes" << endl;
474  }
475 }
476 
477 
479 (
480  const IOobject& io,
481  pointField&& points,
482  faceList&& faces,
483  labelList&& allOwner,
484  labelList&& allNeighbour,
485  const bool syncPar
486 )
487 :
488  polyMesh
489  (
490  io,
491  std::move(points),
492  std::move(faces),
493  std::move(allOwner),
494  std::move(allNeighbour),
495  syncPar
496  ),
497  surfaceInterpolation(*this),
498  data(static_cast<const objectRegistry&>(*this)),
499  boundary_(*this, boundaryMesh()),
500  stitcher_(nullptr),
501  topoChanger_(nullptr),
502  distributor_(nullptr),
503  mover_(nullptr),
504  lduPtr_(nullptr),
505  polyFacesBfPtr_(nullptr),
506  polyBFaceOffsetsPtr_(nullptr),
507  polyBFaceOffsetPatchesPtr_(nullptr),
508  polyBFaceOffsetPatchFacesPtr_(nullptr),
509  polyBFacePatchesPtr_(nullptr),
510  polyBFacePatchFacesPtr_(nullptr),
511  curTimeIndex_(time().timeIndex()),
512  VPtr_(nullptr),
513  V0Ptr_(nullptr),
514  V00Ptr_(nullptr),
515  SfSlicePtr_(nullptr),
516  SfPtr_(nullptr),
517  magSfSlicePtr_(nullptr),
518  magSfPtr_(nullptr),
519  CSlicePtr_(nullptr),
520  CPtr_(nullptr),
521  CfSlicePtr_(nullptr),
522  CfPtr_(nullptr),
523  phiPtr_(nullptr)
524 {
525  if (debug)
526  {
527  Pout<< FUNCTION_NAME << "Constructing fvMesh from components" << endl;
528  }
529 }
530 
531 
533 (
534  const IOobject& io,
535  pointField&& points,
536  faceList&& faces,
537  cellList&& cells,
538  const bool syncPar
539 )
540 :
541  polyMesh
542  (
543  io,
544  std::move(points),
545  std::move(faces),
546  std::move(cells),
547  syncPar
548  ),
549  surfaceInterpolation(*this),
550  data(static_cast<const objectRegistry&>(*this)),
551  boundary_(*this),
552  stitcher_(nullptr),
553  topoChanger_(nullptr),
554  distributor_(nullptr),
555  mover_(nullptr),
556  lduPtr_(nullptr),
557  polyFacesBfPtr_(nullptr),
558  polyBFaceOffsetsPtr_(nullptr),
559  polyBFaceOffsetPatchesPtr_(nullptr),
560  polyBFaceOffsetPatchFacesPtr_(nullptr),
561  polyBFacePatchesPtr_(nullptr),
562  polyBFacePatchFacesPtr_(nullptr),
563  curTimeIndex_(time().timeIndex()),
564  VPtr_(nullptr),
565  V0Ptr_(nullptr),
566  V00Ptr_(nullptr),
567  SfSlicePtr_(nullptr),
568  SfPtr_(nullptr),
569  magSfSlicePtr_(nullptr),
570  magSfPtr_(nullptr),
571  CSlicePtr_(nullptr),
572  CPtr_(nullptr),
573  CfSlicePtr_(nullptr),
574  CfPtr_(nullptr),
575  phiPtr_(nullptr)
576 {
577  if (debug)
578  {
579  Pout<< FUNCTION_NAME << "Constructing fvMesh from components" << endl;
580  }
581 }
582 
583 
585 :
586  polyMesh(Foam::move(mesh)),
587  surfaceInterpolation(Foam::move(mesh)),
588  data(static_cast<const objectRegistry&>(*this)),
589  boundary_(Foam::move(mesh.boundary_)),
590  stitcher_(Foam::move(mesh.stitcher_)),
591  topoChanger_(Foam::move(mesh.topoChanger_)),
592  distributor_(Foam::move(mesh.distributor_)),
593  mover_(Foam::move(mesh.mover_)),
594  lduPtr_(Foam::move(mesh.lduPtr_)),
595  polyFacesBfPtr_(Foam::move(mesh.polyFacesBfPtr_)),
596  polyBFaceOffsetsPtr_(Foam::move(mesh.polyBFaceOffsetsPtr_)),
597  polyBFaceOffsetPatchesPtr_(Foam::move(mesh.polyBFaceOffsetPatchesPtr_)),
598  polyBFaceOffsetPatchFacesPtr_
599  (
600  Foam::move(mesh.polyBFaceOffsetPatchFacesPtr_)
601  ),
602  polyBFacePatchesPtr_(Foam::move(mesh.polyBFacePatchesPtr_)),
603  polyBFacePatchFacesPtr_(Foam::move(mesh.polyBFacePatchFacesPtr_)),
604  curTimeIndex_(mesh.curTimeIndex_),
605  VPtr_(Foam::move(mesh.VPtr_)),
606  V0Ptr_(Foam::move(mesh.V0Ptr_)),
607  V00Ptr_(Foam::move(mesh.V00Ptr_)),
608  SfSlicePtr_(Foam::move(mesh.SfSlicePtr_)),
609  SfPtr_(Foam::move(mesh.SfPtr_)),
610  magSfSlicePtr_(Foam::move(mesh.magSfSlicePtr_)),
611  magSfPtr_(Foam::move(mesh.magSfPtr_)),
612  CSlicePtr_(Foam::move(mesh.CSlicePtr_)),
613  CPtr_(Foam::move(mesh.CPtr_)),
614  CfSlicePtr_(Foam::move(mesh.CfSlicePtr_)),
615  CfPtr_(Foam::move(mesh.CfPtr_)),
616  phiPtr_(Foam::move(mesh.phiPtr_))
617 {
618  if (debug)
619  {
620  Pout<< FUNCTION_NAME << "Moving fvMesh" << endl;
621  }
622 }
623 
624 
625 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
626 
628 {
629  clearOut();
630 }
631 
632 
633 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
634 
636 {
637  return topoChanger_->dynamic();
638 }
639 
640 
642 {
643  return topoChanger_->dynamic() || mover_->dynamic();
644 }
645 
646 
648 {
649  if
650  (
651  stitcher_->stitches()
652  || topoChanger_->dynamic()
653  || distributor_->dynamic()
654  )
655  {
657  }
658 
659  if (!conformal()) stitcher_->disconnect(true, true);
660 
661  const bool hasV00 = V00Ptr_;
662  deleteDemandDrivenData(V00Ptr_);
663 
664  if (!hasV00)
665  {
666  deleteDemandDrivenData(V0Ptr_);
667  }
668 
669  // Set topoChanged_ false before any mesh change
670  topoChanged_ = false;
671  bool updated = topoChanger_->update();
672  topoChanged_ = updated;
673 
674  // Register V0 for distribution
675  if (V0Ptr_)
676  {
677  V0Ptr_->checkIn();
678  }
679 
680  updated = distributor_->update() || updated;
681 
682  // De-register V0 after distribution
683  if (V0Ptr_)
684  {
685  V0Ptr_->checkOut();
686  }
687 
688  if (hasV00)
689  {
690  // If V00 had been set reset to the mapped V0 prior to mesh-motion
691  V00();
692  }
693 
694  return updated;
695 }
696 
697 
699 {
700  if (!conformal()) stitcher_->disconnect(true, true);
701 
702  if (curTimeIndex_ < time().timeIndex() && stitcher_->stitches())
703  {
704  // Store all old-time fields. If we don't do this then we risk
705  // triggering a store in the middle of mapping and potentially
706  // overwriting a mapped old-time field with a not-yet-mapped new-time
707  // field.
708  storeOldTimeFields();
709  }
710 
711  // Do not set moving false
712  // Once the mesh starts moving it is considered to be moving
713  // for the rest of the run
714  const bool moved = mover_->update();
715 
716  curTimeIndex_ = time().timeIndex();
717 
718  stitcher_->connect(true, true, false);
719 
720  return moved;
721 }
722 
723 
725 (
726  const List<polyPatch*> & p,
727  const bool validBoundary
728 )
729 {
730  if (boundary().size())
731  {
733  << " boundary already exists"
734  << abort(FatalError);
735  }
736 
737  // first add polyPatches
738  addPatches(p, validBoundary);
739  boundary_.addPatches(boundaryMesh());
740 }
741 
742 
744 {
745  if (debug)
746  {
747  Pout<< FUNCTION_NAME << "Removing boundary patches." << endl;
748  }
749 
750  // Remove fvBoundaryMesh data first.
751  boundary_.clear();
752  boundary_.setSize(0);
754 
755  clearOut();
756 }
757 
758 
759 void Foam::fvMesh::swap(fvMesh& otherMesh)
760 {
761  // Clear the sliced fields
762  clearGeom();
763 
764  // Clear the current volume and other geometry factors
766 
767  // Clear any non-updateable addressing
768  clearAddressing(true);
769 
770  polyMesh::swap(otherMesh);
771 
772  auto updatePatches = []
773  (
774  const polyPatchList& patches,
775  fvBoundaryMesh& boundaryMesh
776  )
777  {
778  boundaryMesh.setSize(patches.size());
779 
781  {
782  // Construct new processor patches, as the decomposition may have
783  // changed. Leave other patches as is.
784 
785  if (isA<processorPolyPatch>(patches[patchi]))
786  {
787  boundaryMesh.set
788  (
789  patchi,
791  (
792  patches[patchi],
793  boundaryMesh
794  )
795  );
796  }
797  }
798  };
799 
800  updatePatches(boundaryMesh(), boundary_);
801  updatePatches(otherMesh.boundaryMesh(), otherMesh.boundary_);
802 }
803 
804 
806 (
807  const stitchType stitch
808 )
809 {
810  if (debug)
811  {
812  Pout<< FUNCTION_NAME << "Updating fvMesh. ";
813  }
814 
816 
817  if (state == polyMesh::TOPO_PATCH_CHANGE)
818  {
819  boundary_.readUpdate(boundaryMesh());
820  }
821 
822  if
823  (
824  stitcher_.valid()
825  && stitch != stitchType::none
826  && state != polyMesh::UNCHANGED
827  )
828  {
829  stitcher_->disconnect(false, stitch == stitchType::geometric);
830  }
831 
832  if (state == polyMesh::TOPO_PATCH_CHANGE)
833  {
834  if (debug)
835  {
836  Info<< "Boundary and topological update" << endl;
837  }
838 
839  clearOut();
840  }
841  else if (state == polyMesh::TOPO_CHANGE)
842  {
843  if (debug)
844  {
845  Info<< "Topological update" << endl;
846  }
847 
848  clearOut();
849  }
850  else if (state == polyMesh::POINTS_MOVED)
851  {
852  if (debug)
853  {
854  Info<< "Point motion update" << endl;
855  }
856 
857  clearGeom();
858  }
859  else
860  {
861  if (debug)
862  {
863  Info<< "No update" << endl;
864  }
865  }
866 
867  if
868  (
869  stitcher_.valid()
870  && stitch != stitchType::none
871  && state != polyMesh::UNCHANGED
872  )
873  {
874  stitcher_->connect(false, stitch == stitchType::geometric, true);
875  }
876 
877  // If the mesh has been re-stitched with different geometry, then the
878  // finite-volume topology has changed
879  if
880  (
881  stitcher_.valid()
882  && stitcher_->stitches()
883  && state == polyMesh::POINTS_MOVED
884  )
885  {
886  state = polyMesh::TOPO_CHANGE;
887  }
888 
889  return state;
890 }
891 
892 
894 {
895  return boundary_;
896 }
897 
898 
900 {
901  if (!lduPtr_)
902  {
903  lduPtr_ = new fvMeshLduAddressing(*this);
904  }
905 
906  return *lduPtr_;
907 }
908 
909 
911 {
912  return !(polyFacesBfPtr_ && SfPtr_);
913 }
914 
915 
917 {
918  return
919  IOobject
920  (
921  "polyFaces",
922  pointsInstance(),
923  typeName,
924  *this,
925  r,
927  false
928  );
929 }
930 
931 
933 {
934  if (!polyFacesBfPtr_)
935  {
936  polyFacesBfPtr_ =
938  (
939  boundary(),
941  polyFacesPatchTypes(),
942  boundaryMesh().types()
943  );
944  }
945 
946  return *polyFacesBfPtr_;
947 }
948 
949 
952 {
953  if (!polyBFacePatchesPtr_)
954  {
955  const label nPolyBFaces = nFaces() - nInternalFaces();
956 
957  // Count face-poly-bFaces to get the offsets
958  polyBFaceOffsetsPtr_ = new labelList(nPolyBFaces + 1, 0);
959  labelList& offsets = *polyBFaceOffsetsPtr_;
961  {
962  forAll(boundary()[patchi], patchFacei)
963  {
964  const label polyBFacei =
965  (
966  polyFacesBfPtr_
967  ? (*polyFacesBfPtr_)[patchi][patchFacei]
968  : boundary()[patchi].start() + patchFacei
969  )
970  - nInternalFaces();
971 
972  offsets[polyBFacei + 1] ++;
973  }
974  }
975  for (label polyBFacei = 0; polyBFacei < nPolyBFaces; ++ polyBFacei)
976  {
977  offsets[polyBFacei + 1] += offsets[polyBFacei];
978  }
979 
980  // Set the poly-bFace patches and patch-faces, using the offsets as
981  // counters
982  polyBFaceOffsetPatchesPtr_ = new labelList(offsets.last());
983  polyBFaceOffsetPatchFacesPtr_ = new labelList(offsets.last());
984  labelUList& patches = *polyBFaceOffsetPatchesPtr_;
985  labelUList& patchFaces = *polyBFaceOffsetPatchFacesPtr_;
987  {
988  forAll(boundary()[patchi], patchFacei)
989  {
990  const label polyBFacei =
991  (
992  polyFacesBfPtr_
993  ? (*polyFacesBfPtr_)[patchi][patchFacei]
994  : boundary()[patchi].start() + patchFacei
995  )
996  - nInternalFaces();
997  patches[offsets[polyBFacei]] = patchi;
998  patchFaces[offsets[polyBFacei]] = patchFacei;
999  offsets[polyBFacei] ++;
1000  }
1001  }
1002 
1003  // Restore the offsets by removing the count
1004  for
1005  (
1006  label polyBFacei = nPolyBFaces - 1;
1007  polyBFacei >= 0;
1008  -- polyBFacei
1009  )
1010  {
1011  offsets[polyBFacei + 1] = offsets[polyBFacei];
1012  }
1013  offsets[0] = 0;
1014 
1015  // List-lists
1016  polyBFacePatchesPtr_ =
1017  new UCompactListList<label>(offsets, patches);
1018  polyBFacePatchFacesPtr_ =
1019  new UCompactListList<label>(offsets, patchFaces);
1020  }
1021 
1022  return *polyBFacePatchesPtr_;
1023 }
1024 
1025 
1028 {
1029  if (!polyBFacePatchFacesPtr_)
1030  {
1031  polyBFacePatches();
1032  }
1033 
1034  return *polyBFacePatchFacesPtr_;
1035 }
1036 
1037 
1039 {
1040  return stitcher_();
1041 }
1042 
1043 
1045 {
1046  return topoChanger_();
1047 }
1048 
1049 
1051 {
1052  return distributor_();
1053 }
1054 
1055 
1057 {
1058  return mover_();
1059 }
1060 
1061 
1063 {
1064  if (debug)
1065  {
1067  << " nOldCells:" << map.nOldCells()
1068  << " nCells:" << nCells()
1069  << " nOldFaces:" << map.nOldFaces()
1070  << " nFaces:" << nFaces()
1071  << endl;
1072  }
1073 
1074 
1075  // We require geometric properties valid for the old mesh
1076  if
1077  (
1078  map.cellMap().size() != nCells()
1079  || map.faceMap().size() != nFaces()
1080  )
1081  {
1083  << "polyTopoChangeMap does not correspond to the old mesh."
1084  << " nCells:" << nCells()
1085  << " cellMap:" << map.cellMap().size()
1086  << " nOldCells:" << map.nOldCells()
1087  << " nFaces:" << nFaces()
1088  << " faceMap:" << map.faceMap().size()
1089  << " nOldFaces:" << map.nOldFaces()
1090  << exit(FatalError);
1091  }
1092 
1093  // Create a fv mapper
1094  const fvMeshMapper fvMap(*this, map);
1095 
1096  // Map all the volFields in the objectRegistry
1097  #define mapVolFieldType(Type, nullArg) \
1098  MapGeometricFields<Type, fvPatchField, fvMeshMapper, volMesh>(fvMap);
1100 
1101  // Map all the surfaceFields in the objectRegistry
1102  #define mapSurfaceFieldType(Type, nullArg) \
1103  MapGeometricFields<Type, fvsPatchField, fvMeshMapper, surfaceMesh> \
1104  (fvMap);
1106 
1107  // Map all the dimensionedFields in the objectRegistry
1108  #define mapVolInternalFieldType(Type, nullArg) \
1109  MapDimensionedFields<Type, fvMeshMapper, volMesh>(fvMap);
1111 
1112  if (pointMesh::found(*this))
1113  {
1114  // Create the pointMesh mapper
1115  const pointMeshMapper mapper(pointMesh::New(*this), map);
1116 
1117  #define mapPointFieldType(Type, nullArg) \
1118  MapGeometricFields \
1119  < \
1120  Type, \
1121  pointPatchField, \
1122  pointMeshMapper, \
1123  pointMesh \
1124  > \
1125  (mapper);
1127  }
1128 }
1129 
1130 
1132 {
1134 
1135  clearGeom();
1136 
1137  // Update other local data
1138  boundary_.movePoints();
1140 
1141  meshObjects::movePoints<fvMesh>(*this);
1142  meshObjects::movePoints<lduMesh>(*this);
1143 
1144  const_cast<Time&>(time()).functionObjects().movePoints(*this);
1145 }
1146 
1147 
1149 {
1150  // Set moving_ true
1151  // Note: once set it remains true for the rest of the run
1152  moving_ = true;
1153 
1154  // Grab old time volumes if the time has been incremented
1155  // This will update V0, V00
1156  if (curTimeIndex_ < time().timeIndex())
1157  {
1158  storeOldVol(V());
1159  }
1160 
1161  if (!phiPtr_)
1162  {
1163  // Create mesh motion flux
1164  phiPtr_ = new surfaceScalarField
1165  (
1166  IOobject
1167  (
1168  "meshPhi",
1169  this->time().name(),
1170  *this,
1173  true
1174  ),
1175  *this,
1177  );
1178  }
1179  else
1180  {
1181  // Grab old time mesh motion fluxes if the time has been incremented
1182  if (phiPtr_->timeIndex() != time().timeIndex())
1183  {
1184  phiPtr_->oldTime();
1185  }
1186  }
1187 
1188  surfaceScalarField& phi = *phiPtr_;
1189 
1190  // Move the polyMesh and set the mesh motion fluxes to the swept-volumes
1191 
1192  scalar rDeltaT = 1.0/time().deltaTValue();
1193 
1194  tmp<scalarField> tsweptVols = polyMesh::movePoints(p);
1195 
1196  scalarField& sweptVols = tsweptVols.ref();
1197 
1198  phi.primitiveFieldRef() =
1199  scalarField::subField(sweptVols, nInternalFaces());
1200  phi.primitiveFieldRef() *= rDeltaT;
1201 
1202  const fvPatchList& patches = boundary();
1203 
1205 
1207  {
1208  phibf[patchi] = patches[patchi].patchSlice(sweptVols);
1209  phibf[patchi] *= rDeltaT;
1210  }
1211 
1212  // Update or delete the local geometric properties as early as possible so
1213  // they can be used if necessary. These get recreated here instead of
1214  // demand driven since they might do parallel transfers which can conflict
1215  // with when they're actually being used.
1216  // Note that between above "polyMesh::movePoints(p)" and here nothing
1217  // should use the local geometric properties.
1218  updateGeomNotOldVol();
1219 
1220 
1221  // Update other local data
1222  boundary_.movePoints();
1224 
1225  meshObjects::movePoints<fvMesh>(*this);
1226  meshObjects::movePoints<lduMesh>(*this);
1227 
1228  const_cast<Time&>(time()).functionObjects().movePoints(*this);
1229 
1230  return tsweptVols;
1231 }
1232 
1233 
1235 {
1236  // Update polyMesh. This needs to keep volume existent!
1237  polyMesh::topoChange(map);
1238 
1239  if (VPtr_)
1240  {
1241  // Cache old time volumes if they exist and the time has been
1242  // incremented
1243  if (V0Ptr_ && !V0Ptr_->registered())
1244  {
1245  storeOldVol(map.oldCellVolumes());
1246  }
1247 
1248  // Few checks
1249  if (VPtr_ && (VPtr_->size() != map.nOldCells()))
1250  {
1252  << "V:" << V().size()
1253  << " not equal to the number of old cells "
1254  << map.nOldCells()
1255  << exit(FatalError);
1256  }
1257 
1258  if (V0Ptr_ && (V0Ptr_->size() != map.nOldCells()))
1259  {
1261  << "V0:" << V0Ptr_->size()
1262  << " not equal to the number of old cells "
1263  << map.nOldCells()
1264  << exit(FatalError);
1265  }
1266  }
1267 
1268  // Clear the sliced fields
1269  clearGeomNotOldVol();
1270 
1271  // Map the old volume. Just map to new cell labels.
1272  if (V0Ptr_ && !V0Ptr_->registered())
1273  {
1274  const labelList& cellMap = map.cellMap();
1275 
1276  scalarField& V0 = *V0Ptr_;
1277 
1278  scalarField savedV0(V0);
1279  V0.setSize(nCells());
1280 
1281  forAll(V0, i)
1282  {
1283  if (cellMap[i] > -1)
1284  {
1285  V0[i] = savedV0[cellMap[i]];
1286  }
1287  else
1288  {
1289  V0[i] = 0.0;
1290  }
1291  }
1292 
1293  // Inject volume of merged cells
1294  label nMerged = 0;
1295  forAll(map.reverseCellMap(), oldCelli)
1296  {
1297  label index = map.reverseCellMap()[oldCelli];
1298 
1299  if (index < -1)
1300  {
1301  label celli = -index-2;
1302 
1303  V0[celli] += savedV0[oldCelli];
1304 
1305  nMerged++;
1306  }
1307  }
1308 
1309  if (debug)
1310  {
1311  Info<< "Mapping old time volume V0. Merged "
1312  << nMerged << " out of " << nCells() << " cells" << endl;
1313  }
1314  }
1315 
1316  // Map all fields
1317  mapFields(map);
1318 
1319  // Clear the current volume and other geometry factors
1321 
1322  // Clear any non-updateable addressing
1323  clearAddressing(true);
1324 
1325  meshObjects::topoChange<fvMesh>(*this, map);
1326  meshObjects::topoChange<lduMesh>(*this, map);
1327 
1328  const_cast<Time&>(time()).functionObjects().topoChange(map);
1329 
1330  if (topoChanger_.valid())
1331  {
1332  topoChanger_->topoChange(map);
1333  }
1334 
1335  if (distributor_.valid())
1336  {
1337  distributor_->topoChange(map);
1338  }
1339 
1340  if (mover_.valid())
1341  {
1342  mover_->topoChange(map);
1343  }
1344 }
1345 
1346 
1348 {
1349  // Distribute polyMesh data
1350  polyMesh::mapMesh(map);
1351 
1352  // Clear the sliced fields
1353  clearGeomNotOldVol();
1354 
1355  // Clear the current volume and other geometry factors
1357 
1358  // Clear any non-updateable addressing
1359  clearAddressing(true);
1360 
1361  meshObjects::mapMesh<fvMesh>(*this, map);
1362  meshObjects::mapMesh<lduMesh>(*this, map);
1363 
1364  const_cast<Time&>(time()).functionObjects().mapMesh(map);
1365 
1366  topoChanger_->mapMesh(map);
1367  distributor_->mapMesh(map);
1368  mover_->mapMesh(map);
1369 }
1370 
1371 
1373 {
1374  // Distribute polyMesh data
1375  polyMesh::distribute(map);
1376 
1377  // Clear the sliced fields
1378  clearGeomNotOldVol();
1379 
1380  // Clear the current volume and other geometry factors
1382 
1383  // Clear any non-updateable addressing
1384  clearAddressing(true);
1385 
1386  meshObjects::distribute<fvMesh>(*this, map);
1387  meshObjects::distribute<lduMesh>(*this, map);
1388 
1389  const_cast<Time&>(time()).functionObjects().distribute(map);
1390 
1391  topoChanger_->distribute(map);
1392  distributor_->distribute(map);
1393  mover_->distribute(map);
1394 }
1395 
1396 
1398 {
1399  // Clear the geometry fields
1400  clearGeomNotOldVol();
1401 
1402  // Clear the current volume and other geometry factors
1404 
1405  // Clear any non-updateable addressing
1406  clearAddressing(true);
1407 
1408  // Modify the mesh fluxes, if necessary
1409  if (notNull(phi) && phiPtr_)
1410  {
1411  for (label i = 0; i <= phi.nOldTimes(); ++ i)
1412  {
1413  phiRef().oldTime(i) = phi.oldTime(i);
1414  }
1415  }
1416 }
1417 
1418 
1420 (
1421  const surfaceLabelField::Boundary& polyFacesBf,
1422  const surfaceVectorField& Sf,
1423  const surfaceVectorField& Cf,
1424  const surfaceScalarField& phi,
1425  const bool sync
1426 )
1427 {
1428  // Clear the geometry fields
1429  clearGeomNotOldVol();
1430 
1431  // Clear the current volume and other geometry factors
1433 
1434  // Clear any non-updateable addressing
1435  clearAddressing(true);
1436 
1437  // Create non-sliced copies of geometry fields
1438  SfRef();
1439  magSfRef();
1440  CRef();
1441  CfRef();
1442 
1443  // Set the topology
1444  polyFacesBfRef() == polyFacesBf;
1445 
1446  // Set the face geometry
1447  SfRef() == Sf;
1448  magSfRef() == max(mag(Sf), dimensionedScalar(dimArea, rootVSmall));
1449  CRef().boundaryFieldRef() == Cf.boundaryField();
1450  CfRef() == Cf;
1451 
1452  // Communicate processor-coupled cell geometry. Cell-centre processor patch
1453  // fields must contain the (transformed) cell-centre locations on the other
1454  // side of the coupling. This is so that non-conformal patches can
1455  // construct weights and deltas without reference to the poly mesh
1456  // geometry.
1457  //
1458  // Note that the initEvaluate/evaluate communication does a transformation,
1459  // but it is wrong in this instance. A vector field gets transformed as if
1460  // it were a displacement, but the cell-centres need a positional
1461  // transformation. That's why there's the un-transform and re-transform bit
1462  // below just after the evaluate call.
1463  //
1464  // This transform handling is a bit of a hack. It would be nicer to have a
1465  // field attribute which identifies a field as needing a positional
1466  // transformation, and for it to apply automatically within the coupled
1467  // patch field. However, at the moment, the cell centres field is the only
1468  // vol-field containing an absolute position, so the hack is functionally
1469  // sufficient for now.
1470  if (sync && (Pstream::parRun() || !time().processorCase()))
1471  {
1472  volVectorField::Boundary& CBf = CRef().boundaryFieldRef();
1473 
1474  const label nReq = Pstream::nRequests();
1475 
1476  forAll(CBf, patchi)
1477  {
1478  if (isA<processorFvPatch>(CBf[patchi].patch()))
1479  {
1480  CBf[patchi].initEvaluate(Pstream::defaultCommsType);
1481  }
1482  }
1483 
1484  if
1485  (
1486  Pstream::parRun()
1488  )
1489  {
1490  Pstream::waitRequests(nReq);
1491  }
1492 
1493  forAll(CBf, patchi)
1494  {
1495  if (isA<processorFvPatch>(CBf[patchi].patch()))
1496  {
1498 
1499  const transformer& t =
1500  refCast<const processorFvPatch>(CBf[patchi].patch())
1501  .transform();
1502 
1503  t.invTransform(CBf[patchi], CBf[patchi]);
1504  t.transformPosition(CBf[patchi], CBf[patchi]);
1505  }
1506  }
1507  }
1508 
1509  // Modify the mesh fluxes, if necessary
1510  if (notNull(phi) && phiPtr_)
1511  {
1512  for (label i = 0; i <= phi.nOldTimes(); ++ i)
1513  {
1514  phiRef().oldTime(i) = phi.oldTime(i);
1515  }
1516  }
1517 }
1518 
1519 
1521 (
1522  const label insertPatchi,
1523  const polyPatch& patch,
1524  const dictionary& patchFieldDict,
1525  const word& defaultPatchFieldType,
1526  const bool validBoundary
1527 )
1528 {
1529  // Remove my local data (see topoChange)
1530  // Clear mesh motion flux
1531  deleteDemandDrivenData(phiPtr_);
1532 
1533  // Clear the sliced fields
1534  clearGeomNotOldVol();
1535 
1536  // Clear the current volume and other geometry factors
1538 
1539  // Clear any non-updateable addressing
1540  clearAddressing(true);
1541 
1542 
1543  const label sz = boundary_.size();
1544 
1546  (
1547  insertPatchi,
1548  patch,
1549  patchFieldDict,
1550  defaultPatchFieldType,
1551  validBoundary
1552  );
1553 
1554  boundary_.setSize(sz+1);
1555  boundary_.set
1556  (
1557  insertPatchi,
1558  fvPatch::New
1559  (
1560  boundaryMesh()[insertPatchi],
1561  boundary_
1562  )
1563  );
1564 
1565  objectRegistry& db = const_cast<objectRegistry&>(thisDb());
1566  AddPatchFields<volScalarField>
1567  (
1568  db,
1569  insertPatchi,
1570  patchFieldDict,
1571  defaultPatchFieldType,
1572  Zero
1573  );
1574  AddPatchFields<volVectorField>
1575  (
1576  db,
1577  insertPatchi,
1578  patchFieldDict,
1579  defaultPatchFieldType,
1580  Zero
1581  );
1582  AddPatchFields<volSphericalTensorField>
1583  (
1584  db,
1585  insertPatchi,
1586  patchFieldDict,
1587  defaultPatchFieldType,
1588  Zero
1589  );
1590  AddPatchFields<volSymmTensorField>
1591  (
1592  db,
1593  insertPatchi,
1594  patchFieldDict,
1595  defaultPatchFieldType,
1596  Zero
1597  );
1598  AddPatchFields<volTensorField>
1599  (
1600  db,
1601  insertPatchi,
1602  patchFieldDict,
1603  defaultPatchFieldType,
1604  Zero
1605  );
1606 
1607  // Surface fields
1608 
1609  AddPatchFields<surfaceScalarField>
1610  (
1611  db,
1612  insertPatchi,
1613  patchFieldDict,
1614  defaultPatchFieldType,
1615  Zero
1616  );
1617  AddPatchFields<surfaceVectorField>
1618  (
1619  db,
1620  insertPatchi,
1621  patchFieldDict,
1622  defaultPatchFieldType,
1623  Zero
1624  );
1625  AddPatchFields<surfaceSphericalTensorField>
1626  (
1627  db,
1628  insertPatchi,
1629  patchFieldDict,
1630  defaultPatchFieldType,
1631  Zero
1632  );
1633  AddPatchFields<surfaceSymmTensorField>
1634  (
1635  db,
1636  insertPatchi,
1637  patchFieldDict,
1638  defaultPatchFieldType,
1639  Zero
1640  );
1641  AddPatchFields<surfaceTensorField>
1642  (
1643  db,
1644  insertPatchi,
1645  patchFieldDict,
1646  defaultPatchFieldType,
1647  Zero
1648  );
1649 }
1650 
1651 
1653 (
1654  const labelUList& newToOld,
1655  const bool validBoundary
1656 )
1657 {
1658  polyMesh::reorderPatches(newToOld, validBoundary);
1659 
1660  boundary_.shuffle(newToOld, validBoundary);
1661 
1662  objectRegistry& db = const_cast<objectRegistry&>(thisDb());
1663  ReorderPatchFields<volScalarField>(db, newToOld);
1664  ReorderPatchFields<volVectorField>(db, newToOld);
1665  ReorderPatchFields<volSphericalTensorField>(db, newToOld);
1666  ReorderPatchFields<volSymmTensorField>(db, newToOld);
1667  ReorderPatchFields<volTensorField>(db, newToOld);
1668 
1669  ReorderPatchFields<surfaceScalarField>(db, newToOld);
1670  ReorderPatchFields<surfaceVectorField>(db, newToOld);
1671  ReorderPatchFields<surfaceSphericalTensorField>(db, newToOld);
1672  ReorderPatchFields<surfaceSymmTensorField>(db, newToOld);
1673  ReorderPatchFields<surfaceTensorField>(db, newToOld);
1674 }
1675 
1676 
1678 (
1682  const bool write
1683 ) const
1684 {
1685  bool ok = true;
1686 
1687  if (!conformal() && pointsWriteOpt() == IOobject::AUTO_WRITE)
1688  {
1689  // Create a full surface field with the polyFacesBf boundary field to
1690  // write to disk. Make the internal field uniform to save disk space.
1691 
1692  surfaceLabelField polyFaces
1693  (
1694  polyFacesBfIO(IOobject::NO_READ),
1695  *this,
1696  dimless,
1697  labelField(nInternalFaces(), -1),
1698  *polyFacesBfPtr_
1699  );
1700 
1701  ok = ok & polyFaces.write(write);
1702  }
1703 
1704  if (phiPtr_)
1705  {
1706  ok = ok && phiPtr_->write(write);
1707  }
1708 
1709  // Write V0 only if V00 exists
1710  if (V00Ptr_)
1711  {
1712  ok = ok && V0Ptr_->write(write);
1713  }
1714 
1715  if (topoChanger_.valid())
1716  {
1717  topoChanger_->write(write);
1718  }
1719 
1720  if (distributor_.valid())
1721  {
1722  distributor_->write(write);
1723  }
1724 
1725  if (mover_.valid())
1726  {
1727  mover_->write(write);
1728  }
1729 
1730  return ok && polyMesh::writeObject(fmt, ver, cmp, write);
1731 }
1732 
1733 
1734 bool Foam::fvMesh::write(const bool write) const
1735 {
1736  return polyMesh::write(write);
1737 }
1738 
1739 
1740 template<>
1742 Foam::fvMesh::validComponents<Foam::sphericalTensor>() const
1743 {
1745 }
1746 
1747 
1749 {
1750  if (!fvSchemes_.valid())
1751  {
1752  fvSchemes_ = new fvSchemes(*this);
1753  }
1754 
1755  return fvSchemes_;
1756 }
1757 
1758 
1760 {
1761  if (!fvSolution_.valid())
1762  {
1763  fvSolution_ = new fvSolution(*this);
1764  }
1765 
1766  return fvSolution_;
1767 }
1768 
1769 
1770 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
1771 
1772 bool Foam::fvMesh::operator!=(const fvMesh& bm) const
1773 {
1774  return &bm != this;
1775 }
1776 
1777 
1778 bool Foam::fvMesh::operator==(const fvMesh& bm) const
1779 {
1780  return &bm == this;
1781 }
1782 
1783 
1784 // ************************************************************************* //
static const Foam::dimensionedScalar C("C", Foam::dimTemperature, 234.5)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
static bool found(const polyMesh &mesh)
Return true if this DemandDrivenMeshObject is found.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
SubField< scalar > subField
Declare type of subField.
Definition: Field.H:100
static const char *const typeName
Definition: Field.H:105
Generic GeometricBoundaryField class.
void evaluate()
Evaluate boundary conditions.
Generic GeometricField class.
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
GeometricBoundaryField< Type, PatchField, GeoMesh > Boundary
Type of the boundary field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
label nOldTimes() const
Return the number of old time fields stored.
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
readOption
Enumeration defining the read options.
Definition: IOobject.H:117
Version number type.
Definition: IOstream.H:97
streamFormat
Enumeration for the format of data in the stream.
Definition: IOstream.H:87
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:194
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void setSize(const label)
Reset size of List.
Definition: List.C:281
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
Unallocated base class of CompactListList.
T & last()
Return the last element of the list.
Definition: UListI.H:128
static label nRequests()
Get number of outstanding requests.
Definition: UPstream.C:137
static void waitRequests(const label start=0)
Wait until all requests (from start onwards) have finished.
Definition: UPstream.C:147
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
static commsTypes defaultCommsType
Default commsType.
Definition: UPstream.H:272
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
Database for solution and other reduced data.
Definition: data.H:54
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
Foam::fvBoundaryMesh.
Abstract base class for fvMesh movers.
static autoPtr< fvMeshDistributor > New(fvMesh &)
Select, construct and return the fvMeshDistributor.
Foam::fvMeshLduAddressing.
Class holds all the necessary information for mapping fields associated with fvMesh.
Definition: fvMeshMapper.H:56
Abstract base class for fvMesh movers.
Definition: fvMeshMover.H:53
static autoPtr< fvMeshMover > New(fvMesh &)
Select, construct and return the fvMeshMover.
Mesh manipulator that uses the intersection provided by the cyclic non-conformal poly patches to crea...
Abstract base class for fvMesh topology changers.
static autoPtr< fvMeshTopoChanger > New(fvMesh &, const dictionary &dict)
Select, construct and return the fvMeshTopoChanger.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:101
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:406
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:893
const DimensionedField< scalar, volMesh > & V00() const
Return old-old-time cell volumes.
stitchType
Extent to which to stitch on read and readUpdate. By default, a.
Definition: fvMesh.H:113
virtual const lduAddressing & lduAddr() const
Return ldu addressing.
Definition: fvMesh.C:899
virtual void addPatch(const label insertPatchi, const polyPatch &patch, const dictionary &patchFieldDict, const word &defaultPatchFieldType, const bool validBoundary)
Add/insert single patch. If validBoundary the new situation.
Definition: fvMesh.C:1521
const UCompactListList< label > & polyBFacePatches() const
Return poly-bFace-patch addressing.
Definition: fvMesh.C:951
const fvMeshDistributor & distributor() const
Return the distributor function class.
Definition: fvMesh.C:1050
const fvSchemes & schemes() const
Return the fvSchemes.
Definition: fvMesh.C:1748
virtual ~fvMesh()
Destructor.
Definition: fvMesh.C:627
void unconform(const GeometricBoundaryField< label, fvsPatchField, surfaceMesh > &polyFacesBf, const surfaceVectorField &Sf, const surfaceVectorField &Cf, const surfaceScalarField &phi=NullObjectRef< surfaceScalarField >(), const bool sync=true)
Unconform the fvMesh from the polyMesh.
Definition: fvMesh.C:1420
const fvMeshStitcher & stitcher() const
Return the stitcher function class.
Definition: fvMesh.C:1038
const fvSolution & solution() const
Return the fvSchemes.
Definition: fvMesh.C:1759
void addFvPatches(const List< polyPatch * > &, const bool validBoundary=true)
Add boundary patches. Constructor helper.
Definition: fvMesh.C:725
void swap(fvMesh &)
Swap mesh.
Definition: fvMesh.C:759
const GeometricBoundaryField< label, fvsPatchField, surfaceMesh > & polyFacesBf() const
Return face-poly-face addressing.
Definition: fvMesh.C:932
bool operator==(const fvMesh &) const
Definition: fvMesh.C:1778
const UCompactListList< label > & polyBFacePatchFaces() const
Return poly-bFace-patch-face addressing.
Definition: fvMesh.C:1027
bool operator!=(const fvMesh &) const
Definition: fvMesh.C:1772
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool write=true) const
Write the underlying polyMesh and other data.
Definition: fvMesh.C:1678
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Definition: fvMesh.C:1347
void mapFields(const polyTopoChangeMap &map)
Map all fields in time using given map.
Definition: fvMesh.C:1062
virtual bool write(const bool write=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:1734
IOobject polyFacesBfIO(IOobject::readOption r) const
Get the poly faces IO object.
Definition: fvMesh.C:916
bool update()
Update the mesh for topology change, mesh to mesh mapping.
Definition: fvMesh.C:647
virtual void distribute(const polyDistributionMap &map)
Redistribute or update using the given distribution map.
Definition: fvMesh.C:1372
void removeFvBoundary()
Remove boundary patches. Warning: fvPatchFields hold ref to.
Definition: fvMesh.C:743
const word & name() const
Return reference to name.
Definition: fvMesh.H:420
virtual void reorderPatches(const labelUList &newToOld, const bool validBoundary)
Reorder and trim existing patches. If validBoundary the new.
Definition: fvMesh.C:1653
fvMesh(const IOobject &io, const bool changers=true, const stitchType stitch=stitchType::geometric)
Construct from IOobject.
Definition: fvMesh.C:323
bool dynamic() const
Is this mesh dynamic?
Definition: fvMesh.C:641
const fvMeshTopoChanger & topoChanger() const
Return the topo-changer function class.
Definition: fvMesh.C:1044
virtual void topoChange(const polyTopoChangeMap &map)
Update mesh corresponding to the given map.
Definition: fvMesh.C:1234
bool topoChanging() const
Does the mesh topology change?
Definition: fvMesh.C:635
void clearOut()
Clear all geometry and addressing.
Definition: fvMesh.C:274
bool move()
Move the mesh.
Definition: fvMesh.C:698
virtual void setPoints(const pointField &)
Reset the points.
Definition: fvMesh.C:1131
const fvMeshMover & mover() const
Return the mover function class.
Definition: fvMesh.C:1056
void conform(const surfaceScalarField &phi=NullObjectRef< surfaceScalarField >())
Conform the fvMesh to the polyMesh.
Definition: fvMesh.C:1397
bool conformal() const
Return whether the fvMesh is conformal with the polyMesh.
Definition: fvMesh.C:910
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
static autoPtr< fvPatch > New(const polyPatch &, const fvBoundaryMesh &)
Return a pointer to a new patch created on freestore from polyPatch.
Definition: fvPatchNew.C:32
Selector class for finite volume differencing schemes. fvMesh is derived from fvSchemes so that all f...
Definition: fvSchemes.H:53
Selector class for finite volume solution solution. fvMesh is derived from fvSolution so that all fie...
Definition: fvSolution.H:50
The class contains the addressing required by the lduMatrix: upper, lower and losort.
static void clearUpto(objectRegistry &)
Clear all meshObjects derived from FromType up to (but not including)
Registry of regIOobjects.
Traits class for primitives.
Definition: pTraits.H:53
Class holds all the necessary information for mapping fields associated with pointMesh.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
Definition: polyMesh.C:1492
virtual void addPatch(const label insertPatchi, const polyPatch &patch, const dictionary &patchFieldDict, const word &defaultPatchFieldType, const bool validBoundary)
Add/insert single patch. If validBoundary the new situation.
Definition: polyMesh.C:1279
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:91
@ TOPO_PATCH_CHANGE
Definition: polyMesh.H:95
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
void swap(polyMesh &)
Swap mesh.
Definition: polyMesh.C:802
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool write=true) const
Write the underlying polyMesh.
Definition: polyMeshIO.C:519
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
readUpdateState readUpdate()
Update the mesh based on the mesh files saved in.
Definition: polyMeshIO.C:99
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:403
virtual void distribute(const polyDistributionMap &map)
Redistribute or update using the given distribution map.
void removeBoundary()
Remove boundary patches.
Definition: polyMeshClear.C:36
virtual void reorderPatches(const labelUList &newToOld, const bool validBoundary)
Reorder and trim existing patches. If validBoundary the new.
Definition: polyMesh.C:1241
void clearOut()
Clear all geometry and addressing unnecessary for CFD.
virtual void setPoints(const pointField &)
Reset the points.
Definition: polyMesh.C:1449
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
label nOldCells() const
Number of old cells.
const labelList & cellMap() const
Old cell map.
const labelList & reverseCellMap() const
Reverse cell map.
label nOldFaces() const
Number of old faces.
const scalarField & oldCellVolumes() const
const labelList & faceMap() const
Old face map.
void clearAddressing()
Clear topological data.
virtual bool write(const bool write=true) const
Write using setting from DB.
Cell to surface interpolation scheme. Included in fvMesh.
bool movePoints()
Do what is necessary if the mesh has moved.
void clearOut()
Clear all geometry and addressing.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space.
Definition: transformer.H:84
Type invTransform(const Type &) const
Inverse transform the given type.
vector transformPosition(const vector &v) const
Transform the given position.
Definition: transformerI.H:153
Type transform(const Type &) const
Transform the given type.
A class for handling words, derived from string.
Definition: word.H:62
Template functions to aid in the implementation of demand driven data.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
#define nullOldestTimeFields(Type, nullArg)
#define mapVolInternalFieldType(Type, nullArg)
#define mapPointFieldType(Type, nullArg)
#define mapVolFieldType(Type, nullArg)
#define mapSurfaceFieldType(Type, nullArg)
Reconstruct volField from a face flux field.
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField.
label patchi
const pointField & points
const cellShapeList & cells
const fvPatchList & patches
#define FUNCTION_NAME
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
Namespace for OpenFOAM.
bool isFile(const fileName &, const bool checkVariants=true, const bool followLink=true)
Does the name exist as a file in the file system?
Definition: POSIX.C:555
const fileOperation & fileHandler()
Get current file handler.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
static const zero Zero
Definition: zero.H:97
List< label > labelList
A List of labels.
Definition: labelList.H:56
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 & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const dimensionSet dimless
SurfaceField< scalar > surfaceScalarField
messageStream Info
FOR_ALL_FIELD_TYPES(DefineContiguousFvWallLocationDataType)
void deleteDemandDrivenData(DataPtr &dataPtr)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const dimensionSet dimTime
bool notNull(const T &t)
Return true if t is not a reference to the nullObject of type T.
Definition: nullObjectI.H:58
const dimensionSet dimVolume
dimensioned< scalar > mag(const dimensioned< Type > &)
defineTypeNameAndDebug(combustionModel, 0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
Field< label > labelField
Specialisation of Field<T> for label.
Definition: labelField.H:49
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
error FatalError
const dimensionSet dimArea
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label timeIndex
Definition: getTimeIndex.H:4
faceListList boundary(nPatches)
volScalarField & p
Foam::surfaceFields.