All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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-2019 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 "slicedVolFields.H"
30 #include "slicedSurfaceFields.H"
31 #include "SubField.H"
32 #include "demandDrivenData.H"
33 #include "fvMeshLduAddressing.H"
34 #include "mapPolyMesh.H"
35 #include "MapFvFields.H"
36 #include "fvMeshMapper.H"
37 #include "mapClouds.H"
38 #include "MeshObject.H"
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 namespace Foam
43 {
44  defineTypeNameAndDebug(fvMesh, 0);
45 }
46 
47 
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 
50 void Foam::fvMesh::clearGeomNotOldVol()
51 {
52  if (debug)
53  {
54  Pout<< FUNCTION_NAME << "clearGeomNotOldVol" << endl;
55  }
56 
58  <
59  fvMesh,
60  GeometricMeshObject,
61  MoveableMeshObject
62  >(*this);
63 
65  <
66  lduMesh,
67  GeometricMeshObject,
68  MoveableMeshObject
69  >(*this);
70 
72  static_cast<slicedVolScalarField::Internal*>(VPtr_);
74  VPtr_ = nullptr;
75 
76  deleteDemandDrivenData(SfPtr_);
77  deleteDemandDrivenData(magSfPtr_);
79  deleteDemandDrivenData(CfPtr_);
80 }
81 
82 
83 void Foam::fvMesh::updateGeomNotOldVol()
84 {
85  bool haveV = (VPtr_ != nullptr);
86  bool haveSf = (SfPtr_ != nullptr);
87  bool haveMagSf = (magSfPtr_ != nullptr);
88  bool haveCP = (CPtr_ != nullptr);
89  bool haveCf = (CfPtr_ != nullptr);
90 
91  clearGeomNotOldVol();
92 
93  // Now recreate the fields
94  if (haveV)
95  {
96  (void)V();
97  }
98  if (haveSf)
99  {
100  (void)Sf();
101  }
102  if (haveMagSf)
103  {
104  (void)magSf();
105  }
106  if (haveCP)
107  {
108  (void)C();
109  }
110  if (haveCf)
111  {
112  (void)Cf();
113  }
114 }
115 
116 
117 void Foam::fvMesh::clearGeom()
118 {
119  if (debug)
120  {
121  Pout<< FUNCTION_NAME << "Clearing geometric data" << endl;
122  }
123 
124  clearGeomNotOldVol();
125 
126  deleteDemandDrivenData(V0Ptr_);
127  deleteDemandDrivenData(V00Ptr_);
128 
129  // Mesh motion flux cannot be deleted here because the old-time flux
130  // needs to be saved.
131 }
132 
133 
134 void Foam::fvMesh::clearAddressing(const bool isMeshUpdate)
135 {
136  if (debug)
137  {
138  Pout<< FUNCTION_NAME << "isMeshUpdate: " << isMeshUpdate << endl;
139  }
140 
141  if (isMeshUpdate)
142  {
143  // Part of a mesh update. Keep meshObjects that have an updateMesh
144  // callback
146  <
147  fvMesh,
148  TopologicalMeshObject,
149  UpdateableMeshObject
150  >
151  (
152  *this
153  );
155  <
156  lduMesh,
157  TopologicalMeshObject,
158  UpdateableMeshObject
159  >
160  (
161  *this
162  );
163  }
164  else
165  {
166  meshObject::clear<fvMesh, TopologicalMeshObject>(*this);
167  meshObject::clear<lduMesh, TopologicalMeshObject>(*this);
168  }
169  deleteDemandDrivenData(lduPtr_);
170 }
171 
172 
173 void Foam::fvMesh::storeOldVol(const scalarField& V)
174 {
175  if (curTimeIndex_ < time().timeIndex())
176  {
177  if (debug)
178  {
180  << " Storing old time volumes since from time " << curTimeIndex_
181  << " and time now " << time().timeIndex()
182  << " V:" << V.size()
183  << endl;
184  }
185 
186 
187  if (V00Ptr_ && V0Ptr_)
188  {
189  // Copy V0 into V00 storage
190  *V00Ptr_ = *V0Ptr_;
191  }
192 
193  if (V0Ptr_)
194  {
195  // Copy V into V0 storage
196  V0Ptr_->scalarField::operator=(V);
197  }
198  else
199  {
200  // Allocate V0 storage, fill with V
201  V0Ptr_ = new DimensionedField<scalar, volMesh>
202  (
203  IOobject
204  (
205  "V0",
206  time().timeName(),
207  *this,
210  false
211  ),
212  *this,
213  dimVolume
214  );
215  scalarField& V0 = *V0Ptr_;
216  // Note: V0 now sized with current mesh, not with (potentially
217  // different size) V.
218  V0.setSize(V.size());
219  V0 = V;
220  }
221 
222  curTimeIndex_ = time().timeIndex();
223 
224  if (debug)
225  {
227  << " Stored old time volumes V0:" << V0Ptr_->size()
228  << endl;
229  if (V00Ptr_)
230  {
232  << " Stored oldold time volumes V00:" << V00Ptr_->size()
233  << endl;
234  }
235  }
236  }
237 }
238 
239 
241 {
242  clearGeom();
244 
245  clearAddressing();
246 
247  // Clear mesh motion flux
248  deleteDemandDrivenData(phiPtr_);
249 
251 }
252 
253 
254 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
255 
257 :
258  polyMesh(io),
259  surfaceInterpolation(*this),
260  fvSchemes(static_cast<const objectRegistry&>(*this)),
261  fvSolution(static_cast<const objectRegistry&>(*this)),
262  data(static_cast<const objectRegistry&>(*this)),
263  boundary_(*this, boundaryMesh()),
264  lduPtr_(nullptr),
265  curTimeIndex_(time().timeIndex()),
266  VPtr_(nullptr),
267  V0Ptr_(nullptr),
268  V00Ptr_(nullptr),
269  SfPtr_(nullptr),
270  magSfPtr_(nullptr),
271  CPtr_(nullptr),
272  CfPtr_(nullptr),
273  phiPtr_(nullptr)
274 {
275  if (debug)
276  {
277  Pout<< FUNCTION_NAME << "Constructing fvMesh from IOobject" << endl;
278  }
279 
280  // Check the existence of the cell volumes and read if present
281  // and set the storage of V00
282  if (fileHandler().isFile(time().timePath()/"V0"))
283  {
285  (
286  IOobject
287  (
288  "V0",
289  time().timeName(),
290  *this,
293  false
294  ),
295  *this
296  );
297 
298  V00();
299  }
300 
301  // Check the existence of the mesh fluxes and read if present
302  if (fileHandler().isFile(time().timePath()/"meshPhi"))
303  {
304  phiPtr_ = new surfaceScalarField
305  (
306  IOobject
307  (
308  "meshPhi",
309  time().timeName(),
310  *this,
313  true
314  ),
315  *this
316  );
317  }
318 }
319 
320 
322 (
323  const IOobject& io,
324  pointField&& points,
325  const cellShapeList& shapes,
326  const faceListList& boundaryFaces,
327  const wordList& boundaryPatchNames,
328  const PtrList<dictionary>& boundaryDicts,
329  const word& defaultBoundaryPatchName,
330  const word& defaultBoundaryPatchType,
331  const bool syncPar
332 )
333 :
334  polyMesh
335  (
336  io,
337  move(points),
338  shapes,
339  boundaryFaces,
340  boundaryPatchNames,
341  boundaryDicts,
342  defaultBoundaryPatchName,
343  defaultBoundaryPatchType,
344  syncPar
345  ),
346  surfaceInterpolation(*this),
347  fvSchemes(static_cast<const objectRegistry&>(*this)),
348  fvSolution(static_cast<const objectRegistry&>(*this)),
349  data(static_cast<const objectRegistry&>(*this)),
350  boundary_(*this, boundaryMesh()),
351  lduPtr_(nullptr),
352  curTimeIndex_(time().timeIndex()),
353  VPtr_(nullptr),
354  V0Ptr_(nullptr),
355  V00Ptr_(nullptr),
356  SfPtr_(nullptr),
357  magSfPtr_(nullptr),
358  CPtr_(nullptr),
359  CfPtr_(nullptr),
360  phiPtr_(nullptr)
361 {
362  if (debug)
363  {
364  Pout<< FUNCTION_NAME << "Constructing fvMesh from cellShapes" << endl;
365  }
366 }
367 
368 
370 (
371  const IOobject& io,
372  pointField&& points,
373  faceList&& faces,
374  labelList&& allOwner,
375  labelList&& allNeighbour,
376  const bool syncPar
377 )
378 :
379  polyMesh
380  (
381  io,
382  move(points),
383  move(faces),
384  move(allOwner),
385  move(allNeighbour),
386  syncPar
387  ),
388  surfaceInterpolation(*this),
389  fvSchemes(static_cast<const objectRegistry&>(*this)),
390  fvSolution(static_cast<const objectRegistry&>(*this)),
391  data(static_cast<const objectRegistry&>(*this)),
392  boundary_(*this, boundaryMesh()),
393  lduPtr_(nullptr),
394  curTimeIndex_(time().timeIndex()),
395  VPtr_(nullptr),
396  V0Ptr_(nullptr),
397  V00Ptr_(nullptr),
398  SfPtr_(nullptr),
399  magSfPtr_(nullptr),
400  CPtr_(nullptr),
401  CfPtr_(nullptr),
402  phiPtr_(nullptr)
403 {
404  if (debug)
405  {
406  Pout<< FUNCTION_NAME << "Constructing fvMesh from components" << endl;
407  }
408 }
409 
410 
412 (
413  const IOobject& io,
414  pointField&& points,
415  faceList&& faces,
416  cellList&& cells,
417  const bool syncPar
418 )
419 :
420  polyMesh(io, move(points), move(faces), move(cells), syncPar),
421  surfaceInterpolation(*this),
422  fvSchemes(static_cast<const objectRegistry&>(*this)),
423  fvSolution(static_cast<const objectRegistry&>(*this)),
424  data(static_cast<const objectRegistry&>(*this)),
425  boundary_(*this),
426  lduPtr_(nullptr),
427  curTimeIndex_(time().timeIndex()),
428  VPtr_(nullptr),
429  V0Ptr_(nullptr),
430  V00Ptr_(nullptr),
431  SfPtr_(nullptr),
432  magSfPtr_(nullptr),
433  CPtr_(nullptr),
434  CfPtr_(nullptr),
435  phiPtr_(nullptr)
436 {
437  if (debug)
438  {
439  Pout<< FUNCTION_NAME << "Constructing fvMesh from components" << endl;
440  }
441 }
442 
443 
444 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
445 
447 {
448  clearOut();
449 }
450 
451 
452 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
453 
455 (
456  const List<polyPatch*> & p,
457  const bool validBoundary
458 )
459 {
460  if (boundary().size())
461  {
463  << " boundary already exists"
464  << abort(FatalError);
465  }
466 
467  // first add polyPatches
468  addPatches(p, validBoundary);
469  boundary_.addPatches(boundaryMesh());
470 }
471 
472 
474 {
475  if (debug)
476  {
477  Pout<< FUNCTION_NAME << "Removing boundary patches." << endl;
478  }
479 
480  // Remove fvBoundaryMesh data first.
481  boundary_.clear();
482  boundary_.setSize(0);
484 
485  clearOut();
486 }
487 
488 
490 {
491  if (debug)
492  {
493  Pout<< FUNCTION_NAME << "Updating fvMesh. ";
494  }
495 
497 
498  if (state == polyMesh::TOPO_PATCH_CHANGE)
499  {
500  if (debug)
501  {
502  Info<< "Boundary and topological update" << endl;
503  }
504 
505  boundary_.readUpdate(boundaryMesh());
506 
507  clearOut();
508 
509  }
510  else if (state == polyMesh::TOPO_CHANGE)
511  {
512  if (debug)
513  {
514  Info<< "Topological update" << endl;
515  }
516 
517  clearOut();
518  }
519  else if (state == polyMesh::POINTS_MOVED)
520  {
521  if (debug)
522  {
523  Info<< "Point motion update" << endl;
524  }
525 
526  clearGeom();
527  }
528  else
529  {
530  if (debug)
531  {
532  Info<< "No update" << endl;
533  }
534  }
535 
536  return state;
537 }
538 
539 
541 {
542  return boundary_;
543 }
544 
545 
547 {
548  if (!lduPtr_)
549  {
550  lduPtr_ = new fvMeshLduAddressing(*this);
551  }
552 
553  return *lduPtr_;
554 }
555 
556 
558 {
559  if (debug)
560  {
562  << " nOldCells:" << meshMap.nOldCells()
563  << " nCells:" << nCells()
564  << " nOldFaces:" << meshMap.nOldFaces()
565  << " nFaces:" << nFaces()
566  << endl;
567  }
568 
569 
570  // We require geometric properties valid for the old mesh
571  if
572  (
573  meshMap.cellMap().size() != nCells()
574  || meshMap.faceMap().size() != nFaces()
575  )
576  {
578  << "mapPolyMesh does not correspond to the old mesh."
579  << " nCells:" << nCells()
580  << " cellMap:" << meshMap.cellMap().size()
581  << " nOldCells:" << meshMap.nOldCells()
582  << " nFaces:" << nFaces()
583  << " faceMap:" << meshMap.faceMap().size()
584  << " nOldFaces:" << meshMap.nOldFaces()
585  << exit(FatalError);
586  }
587 
588  // Create a mapper
589  const fvMeshMapper mapper(*this, meshMap);
590 
591  // Map all the volFields in the objectRegistry
592  MapGeometricFields<scalar, fvPatchField, fvMeshMapper, volMesh>
593  (mapper);
594  MapGeometricFields<vector, fvPatchField, fvMeshMapper, volMesh>
595  (mapper);
596  MapGeometricFields<sphericalTensor, fvPatchField, fvMeshMapper, volMesh>
597  (mapper);
598  MapGeometricFields<symmTensor, fvPatchField, fvMeshMapper, volMesh>
599  (mapper);
600  MapGeometricFields<tensor, fvPatchField, fvMeshMapper, volMesh>
601  (mapper);
602 
603  // Map all the surfaceFields in the objectRegistry
604  MapGeometricFields<scalar, fvsPatchField, fvMeshMapper, surfaceMesh>
605  (mapper);
606  MapGeometricFields<vector, fvsPatchField, fvMeshMapper, surfaceMesh>
607  (mapper);
608  MapGeometricFields<symmTensor, fvsPatchField, fvMeshMapper, surfaceMesh>
609  (mapper);
610  MapGeometricFields<symmTensor, fvsPatchField, fvMeshMapper, surfaceMesh>
611  (mapper);
612  MapGeometricFields<tensor, fvsPatchField, fvMeshMapper, surfaceMesh>
613  (mapper);
614 
615  // Map all the dimensionedFields in the objectRegistry
616  MapDimensionedFields<scalar, fvMeshMapper, volMesh>(mapper);
617  MapDimensionedFields<vector, fvMeshMapper, volMesh>(mapper);
618  MapDimensionedFields<sphericalTensor, fvMeshMapper, volMesh>(mapper);
619  MapDimensionedFields<symmTensor, fvMeshMapper, volMesh>(mapper);
620  MapDimensionedFields<tensor, fvMeshMapper, volMesh>(mapper);
621 
622  // Map all the clouds in the objectRegistry
623  mapClouds(*this, meshMap);
624 
625 
626  const labelList& cellMap = meshMap.cellMap();
627 
628  // Map the old volume. Just map to new cell labels.
629  if (V0Ptr_)
630  {
631  scalarField& V0 = *V0Ptr_;
632 
633  scalarField savedV0(V0);
634  V0.setSize(nCells());
635 
636  forAll(V0, i)
637  {
638  if (cellMap[i] > -1)
639  {
640  V0[i] = savedV0[cellMap[i]];
641  }
642  else
643  {
644  V0[i] = 0.0;
645  }
646  }
647 
648  // Inject volume of merged cells
649  label nMerged = 0;
650  forAll(meshMap.reverseCellMap(), oldCelli)
651  {
652  label index = meshMap.reverseCellMap()[oldCelli];
653 
654  if (index < -1)
655  {
656  label celli = -index-2;
657 
658  V0[celli] += savedV0[oldCelli];
659 
660  nMerged++;
661  }
662  }
663 
664  if (debug)
665  {
666  Info<< "Mapping old time volume V0. Merged "
667  << nMerged << " out of " << nCells() << " cells" << endl;
668  }
669  }
670 
671 
672  // Map the old-old volume. Just map to new cell labels.
673  if (V00Ptr_)
674  {
675  scalarField& V00 = *V00Ptr_;
676 
677  scalarField savedV00(V00);
678  V00.setSize(nCells());
679 
680  forAll(V00, i)
681  {
682  if (cellMap[i] > -1)
683  {
684  V00[i] = savedV00[cellMap[i]];
685  }
686  else
687  {
688  V00[i] = 0.0;
689  }
690  }
691 
692  // Inject volume of merged cells
693  label nMerged = 0;
694  forAll(meshMap.reverseCellMap(), oldCelli)
695  {
696  label index = meshMap.reverseCellMap()[oldCelli];
697 
698  if (index < -1)
699  {
700  label celli = -index-2;
701 
702  V00[celli] += savedV00[oldCelli];
703  nMerged++;
704  }
705  }
706 
707  if (debug)
708  {
709  Info<< "Mapping old time volume V00. Merged "
710  << nMerged << " out of " << nCells() << " cells" << endl;
711  }
712  }
713 }
714 
715 
717 {
718  // Grab old time volumes if the time has been incremented
719  // This will update V0, V00
720  if (curTimeIndex_ < time().timeIndex())
721  {
722  storeOldVol(V());
723  }
724 
725  if (!phiPtr_)
726  {
727  // Create mesh motion flux
728  phiPtr_ = new surfaceScalarField
729  (
730  IOobject
731  (
732  "meshPhi",
733  this->time().timeName(),
734  *this,
737  true
738  ),
739  *this,
741  );
742  }
743  else
744  {
745  // Grab old time mesh motion fluxes if the time has been incremented
746  if (phiPtr_->timeIndex() != time().timeIndex())
747  {
748  phiPtr_->oldTime();
749  }
750  }
751 
752  surfaceScalarField& phi = *phiPtr_;
753 
754  // Move the polyMesh and set the mesh motion fluxes to the swept-volumes
755 
756  scalar rDeltaT = 1.0/time().deltaTValue();
757 
758  tmp<scalarField> tsweptVols = polyMesh::movePoints(p);
759  scalarField& sweptVols = tsweptVols.ref();
760 
761  phi.primitiveFieldRef() =
763  phi.primitiveFieldRef() *= rDeltaT;
764 
765  const fvPatchList& patches = boundary();
766 
767  surfaceScalarField::Boundary& phibf = phi.boundaryFieldRef();
768 
769  forAll(patches, patchi)
770  {
771  phibf[patchi] = patches[patchi].patchSlice(sweptVols);
772  phibf[patchi] *= rDeltaT;
773  }
774 
775  // Update or delete the local geometric properties as early as possible so
776  // they can be used if necessary. These get recreated here instead of
777  // demand driven since they might do parallel transfers which can conflict
778  // with when they're actually being used.
779  // Note that between above "polyMesh::movePoints(p)" and here nothing
780  // should use the local geometric properties.
781  updateGeomNotOldVol();
782 
783 
784  // Update other local data
785  boundary_.movePoints();
787 
788  meshObject::movePoints<fvMesh>(*this);
789  meshObject::movePoints<lduMesh>(*this);
790 
791  return tsweptVols;
792 }
793 
794 
796 {
797  // Update polyMesh. This needs to keep volume existent!
799 
800  if (VPtr_)
801  {
802  // Grab old time volumes if the time has been incremented
803  // This will update V0, V00
804  storeOldVol(mpm.oldCellVolumes());
805 
806  // Few checks
807  if (VPtr_ && (V().size() != mpm.nOldCells()))
808  {
810  << "V:" << V().size()
811  << " not equal to the number of old cells "
812  << mpm.nOldCells()
813  << exit(FatalError);
814  }
815  if (V0Ptr_ && (V0Ptr_->size() != mpm.nOldCells()))
816  {
818  << "V0:" << V0Ptr_->size()
819  << " not equal to the number of old cells "
820  << mpm.nOldCells()
821  << exit(FatalError);
822  }
823  if (V00Ptr_ && (V00Ptr_->size() != mpm.nOldCells()))
824  {
826  << "V0:" << V00Ptr_->size()
827  << " not equal to the number of old cells "
828  << mpm.nOldCells()
829  << exit(FatalError);
830  }
831  }
832 
833 
834  // Clear mesh motion flux (note: could instead save & map like volumes)
835  // deleteDemandDrivenData(phiPtr_);
836 
837  // Clear the sliced fields
838  clearGeomNotOldVol();
839 
840  // Map all fields
841  mapFields(mpm);
842 
843  // Clear the current volume and other geometry factors
845 
846  // Clear any non-updateable addressing
847  clearAddressing(true);
848 
849  meshObject::updateMesh<fvMesh>(*this, mpm);
850  meshObject::updateMesh<lduMesh>(*this, mpm);
851 }
852 
853 
855 (
856  const label insertPatchi,
857  const polyPatch& patch,
858  const dictionary& patchFieldDict,
859  const word& defaultPatchFieldType,
860  const bool validBoundary
861 )
862 {
863  // Remove my local data (see updateMesh)
864  // Clear mesh motion flux
865  deleteDemandDrivenData(phiPtr_);
866 
867  // Clear the sliced fields
868  clearGeomNotOldVol();
869 
870  // Clear the current volume and other geometry factors
872 
873  // Clear any non-updateable addressing
874  clearAddressing(true);
875 
876 
877  const label sz = boundary_.size();
878 
880  (
881  insertPatchi,
882  patch,
883  patchFieldDict,
884  defaultPatchFieldType,
885  validBoundary
886  );
887 
888  boundary_.setSize(sz+1);
889  boundary_.set
890  (
891  insertPatchi,
893  (
894  boundaryMesh()[insertPatchi],
895  boundary_
896  )
897  );
898 
899  objectRegistry& db = const_cast<objectRegistry&>(thisDb());
900  AddPatchFields<volScalarField>
901  (
902  db,
903  insertPatchi,
904  patchFieldDict,
905  defaultPatchFieldType,
906  Zero
907  );
908  AddPatchFields<volVectorField>
909  (
910  db,
911  insertPatchi,
912  patchFieldDict,
913  defaultPatchFieldType,
914  Zero
915  );
916  AddPatchFields<volSphericalTensorField>
917  (
918  db,
919  insertPatchi,
920  patchFieldDict,
921  defaultPatchFieldType,
922  Zero
923  );
924  AddPatchFields<volSymmTensorField>
925  (
926  db,
927  insertPatchi,
928  patchFieldDict,
929  defaultPatchFieldType,
930  Zero
931  );
932  AddPatchFields<volTensorField>
933  (
934  db,
935  insertPatchi,
936  patchFieldDict,
937  defaultPatchFieldType,
938  Zero
939  );
940 
941  // Surface fields
942 
943  AddPatchFields<surfaceScalarField>
944  (
945  db,
946  insertPatchi,
947  patchFieldDict,
948  defaultPatchFieldType,
949  Zero
950  );
951  AddPatchFields<surfaceVectorField>
952  (
953  db,
954  insertPatchi,
955  patchFieldDict,
956  defaultPatchFieldType,
957  Zero
958  );
959  AddPatchFields<surfaceSphericalTensorField>
960  (
961  db,
962  insertPatchi,
963  patchFieldDict,
964  defaultPatchFieldType,
965  Zero
966  );
967  AddPatchFields<surfaceSymmTensorField>
968  (
969  db,
970  insertPatchi,
971  patchFieldDict,
972  defaultPatchFieldType,
973  Zero
974  );
975  AddPatchFields<surfaceTensorField>
976  (
977  db,
978  insertPatchi,
979  patchFieldDict,
980  defaultPatchFieldType,
981  Zero
982  );
983 }
984 
985 
987 (
988  const labelUList& newToOld,
989  const bool validBoundary
990 )
991 {
992  polyMesh::reorderPatches(newToOld, validBoundary);
993 
994  boundary_.shuffle(newToOld, validBoundary);
995 
996  objectRegistry& db = const_cast<objectRegistry&>(thisDb());
997  ReorderPatchFields<volScalarField>(db, newToOld);
998  ReorderPatchFields<volVectorField>(db, newToOld);
999  ReorderPatchFields<volSphericalTensorField>(db, newToOld);
1000  ReorderPatchFields<volSymmTensorField>(db, newToOld);
1001  ReorderPatchFields<volTensorField>(db, newToOld);
1002 
1003  ReorderPatchFields<surfaceScalarField>(db, newToOld);
1004  ReorderPatchFields<surfaceVectorField>(db, newToOld);
1005  ReorderPatchFields<surfaceSphericalTensorField>(db, newToOld);
1006  ReorderPatchFields<surfaceSymmTensorField>(db, newToOld);
1007  ReorderPatchFields<surfaceTensorField>(db, newToOld);
1008 }
1009 
1010 
1016  const bool write
1017 ) const
1018 {
1019  bool ok = true;
1020  if (phiPtr_)
1021  {
1022  ok = phiPtr_->write(write);
1023  }
1024 
1025  // Write V0 only if V00 exists
1026  if (V00Ptr_)
1027  {
1028  ok = ok && V0Ptr_->write(write);
1029  }
1030 
1031  return ok && polyMesh::writeObject(fmt, ver, cmp, write);
1032 }
1033 
1034 
1035 bool Foam::fvMesh::write(const bool write) const
1036 {
1037  return polyMesh::write(write);
1038 }
1039 
1040 
1041 template<>
1043 Foam::fvMesh::validComponents<Foam::sphericalTensor>() const
1044 {
1046 }
1047 
1048 
1049 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
1050 
1051 bool Foam::fvMesh::operator!=(const fvMesh& bm) const
1052 {
1053  return &bm != this;
1054 }
1055 
1056 
1057 bool Foam::fvMesh::operator==(const fvMesh& bm) const
1058 {
1059  return &bm == this;
1060 }
1061 
1062 
1063 // ************************************************************************* //
Foam::surfaceFields.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
virtual void reorderPatches(const labelUList &newToOld, const bool validBoundary)
Reorder and trim existing patches. If validBoundary the new.
Definition: polyMesh.C:1011
void clearOut()
Clear all geometry and addressing.
Definition: fvMesh.C:240
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
Definition: polyMesh.C:1287
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
void clearAddressing()
Clear topological data.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
void removeFvBoundary()
Remove boundary patches. Warning: fvPatchFields hold ref to.
Definition: fvMesh.C:473
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
label nOldCells() const
Number of old cells.
Definition: mapPolyMesh.H:376
bool set(const label) const
Is element set.
Definition: PtrListI.H:65
const labelList & cellMap() const
Old cell map.
Definition: mapPolyMesh.H:424
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void mapClouds(const objectRegistry &db, const mapPolyMesh &mapper)
Generic Geometric field mapper.
Definition: mapClouds.H:48
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
fvSchemes(const objectRegistry &obr)
Construct for objectRegistry.
Definition: fvSchemes.C:195
void addFvPatches(const List< polyPatch *> &, const bool validBoundary=true)
Add boundary patches. Constructor helper.
Definition: fvMesh.C:455
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static int debug
Debug switch.
Definition: fvSchemes.H:97
label nInternalFaces() const
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:855
label nFaces() const
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
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
void clearOut()
Clear all geometry and addressing unnecessary for CFD.
virtual void reorderPatches(const labelUList &newToOld, const bool validBoundary)
Reorder and trim existing patches. If validBoundary the new.
Definition: fvMesh.C:987
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
IOobject(const word &name, const fileName &instance, const objectRegistry &registry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Definition: IOobject.C:177
label nCells() const
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
const surfaceScalarField & phi() const
Return cell face motion fluxes.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Traits class for primitives.
Definition: pTraits.H:50
bool movePoints()
Do what is necessary if the mesh has moved.
const DimensionedField< scalar, volMesh > & V00() const
Return old-old-time cell volumes.
const cellList & cells() const
Cell to surface interpolation scheme. Included in fvMesh.
patches[0]
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
virtual bool write(const bool write=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:1035
DimensionedField< Type, GeoMesh > Internal
Type of the internal field from which this GeometricField is derived.
Class holds all the necessary information for mapping fields associated with fvMesh.
Definition: fvMeshMapper.H:55
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
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:1012
virtual const objectRegistry & thisDb() const
Return the object registry - resolve conflict polyMesh/lduMesh.
Definition: fvMesh.H:245
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
virtual ~fvMesh()
Destructor.
Definition: fvMesh.C:446
void readUpdate(const polyBoundaryMesh &)
Update boundary based on new polyBoundaryMesh.
bool operator!=(const fvMesh &) const
Definition: fvMesh.C:1051
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
fvSolution(const objectRegistry &obr)
Construct for objectRegistry.
Definition: fvSolution.H:56
virtual void updateMesh(const mapPolyMesh &mpm)
Update mesh corresponding to the given map.
Definition: fvMesh.C:795
virtual readUpdateState readUpdate()
Update the mesh based on the mesh files saved in time.
Definition: fvMesh.C:489
virtual void updateMesh(const mapPolyMesh &mpm)
Update the mesh corresponding to given map.
surfaceInterpolation(const fvMesh &)
Construct given an fvMesh.
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual const lduAddressing & lduAddr() const
Return ldu addressing.
Definition: fvMesh.C:546
static void clearUpto(objectRegistry &)
Clear all meshObject derived from FromType up to (but not including)
Definition: MeshObject.C:494
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:1049
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
streamFormat
Enumeration for the format of data in the stream.
Definition: IOstream.H:86
void shuffle(const labelUList &newToOld, const bool validBoundary)
Reorders patches. Ordering does not have to be done in.
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:97
const fileOperation & fileHandler()
Get current file handler.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
errorManip< error > abort(error &err)
Definition: errorManip.H:131
void addPatches(const List< polyPatch *> &, const bool validBoundary=true)
Add boundary patches.
Definition: polyMesh.C:909
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:193
void removeBoundary()
Remove boundary patches.
Definition: polyMeshClear.C:36
Foam::fvMeshLduAddressing.
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
const labelList & reverseCellMap() const
Reverse cell map.
Definition: mapPolyMesh.H:521
#define FUNCTION_NAME
defineTypeNameAndDebug(combustionModel, 0)
Database for solution and other reduced data.
Definition: data.H:51
const labelList & faceMap() const
Old face map.
Definition: mapPolyMesh.H:399
Generic Geometric field mapper. For "real" mapping, add template specialisations for mapping of inter...
label timeIndex() const
Return the time index of the field.
fvMesh(const IOobject &io)
Construct from IOobject.
Definition: fvMesh.C:256
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
void setSize(const label)
Reset size of List.
Definition: List.C:281
Template functions to aid in the implementation of demand driven data.
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
label patchi
const scalarField & oldCellVolumes() const
Definition: mapPolyMesh.H:640
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
Selector class for finite volume solution solution. fvMesh is derived from fvSolution so that all fie...
Definition: fvSolution.H:47
void clearOut()
Clear all geometry and addressing.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
bool operator==(const fvMesh &) const
Definition: fvMesh.C:1057
Foam::fvBoundaryMesh.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
Selector class for finite volume differencing schemes. fvMesh is derived from fvSchemes so that all f...
Definition: fvSchemes.H:50
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label timeIndex
Definition: getTimeIndex.H:4
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
Version number type.
Definition: IOstream.H:96
void movePoints()
Correct patches after moving points.
void clear()
Clear the PtrList, i.e. set size to zero deleting all the.
Definition: PtrList.C:174
messageStream Info
SubField< scalar > subField
Declare type of subField.
Definition: Field.H:100
The class contains the addressing required by the lduMatrix: upper, lower and losort.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
static autoPtr< fvPatch > New(const polyPatch &, const fvBoundaryMesh &)
Return a pointer to a new patch created on freestore from polyPatch.
Definition: fvPatchNew.C:32
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:88
virtual bool write(const bool write=true) const
Write using setting from DB.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual readUpdateState readUpdate()
Update the mesh based on the mesh files saved in.
Definition: polyMeshIO.C:71
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Registry of regIOobjects.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:322
void deleteDemandDrivenData(DataPtr &dataPtr)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool write) const
Write the objects.
label nOldFaces() const
Number of old faces.
Definition: mapPolyMesh.H:370
Namespace for OpenFOAM.
virtual void mapFields(const mapPolyMesh &mpm)
Map all fields in time using given map.
Definition: fvMesh.C:557
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540
polyMesh(const IOobject &io)
Construct from IOobject.
Definition: polyMesh.C:163