Go to the documentation of this file.
54 void Foam::fvMesh::clearGeomNotOldVol()
70 slicedVolScalarField::Internal* VPtr =
71 static_cast<slicedVolScalarField::Internal*
>(VPtr_);
82 void Foam::fvMesh::updateGeomNotOldVol()
84 bool haveV = (VPtr_ !=
nullptr);
85 bool haveSf = (SfPtr_ !=
nullptr);
86 bool haveMagSf = (magSfPtr_ !=
nullptr);
87 bool haveCP = (CPtr_ !=
nullptr);
88 bool haveCf = (CfPtr_ !=
nullptr);
116 void Foam::fvMesh::clearGeom()
118 clearGeomNotOldVol();
142 TopologicalMeshObject,
151 TopologicalMeshObject,
160 meshObject::clear<fvMesh, TopologicalMeshObject>(*
this);
161 meshObject::clear<lduMesh, TopologicalMeshObject>(*
this);
167 void Foam::fvMesh::storeOldVol(
const scalarField& V)
174 <<
" Storing old time volumes since from time " << curTimeIndex_
175 <<
" and time now " << time().timeIndex()
181 if (V00Ptr_ && V0Ptr_)
190 V0Ptr_->scalarField::operator=(V);
195 V0Ptr_ =
new DimensionedField<scalar, volMesh>
212 V0.setSize(V.size());
216 curTimeIndex_ = time().timeIndex();
221 <<
" Stored old time volumes V0:" << V0Ptr_->size()
226 <<
" Stored oldold time volumes V00:" << V00Ptr_->size()
354 std::move(allNeighbour),
420 Foam::fvMesh::fvMesh(
const IOobject& io,
const zero,
const bool syncPar)
494 const bool validBoundary
500 <<
" boundary already exists"
504 addPatches(plist, validBoundary);
512 const bool validBoundary
518 addFvPatches(plist, validBoundary);
531 boundary_.setSize(0);
551 Info<<
"Boundary and topological update" <<
endl;
563 Info<<
"Topological update" <<
endl;
572 Info<<
"Point motion update" <<
endl;
602 <<
" calculating fvMeshLduAddressing from nFaces:"
619 <<
" nCells:" << nCells()
621 <<
" nFaces:" << nFaces()
629 meshMap.
cellMap().size() != nCells()
630 || meshMap.
faceMap().size() != nFaces()
634 <<
"mapPolyMesh does not correspond to the old mesh."
635 <<
" nCells:" << nCells()
636 <<
" cellMap:" << meshMap.
cellMap().size()
638 <<
" nFaces:" << nFaces()
639 <<
" faceMap:" << meshMap.
faceMap().size()
648 MapGeometricFields<scalar, fvPatchField, fvMeshMapper, volMesh>
650 MapGeometricFields<vector, fvPatchField, fvMeshMapper, volMesh>
652 MapGeometricFields<sphericalTensor, fvPatchField, fvMeshMapper, volMesh>
654 MapGeometricFields<symmTensor, fvPatchField, fvMeshMapper, volMesh>
656 MapGeometricFields<tensor, fvPatchField, fvMeshMapper, volMesh>
660 MapGeometricFields<scalar, fvsPatchField, fvMeshMapper, surfaceMesh>
662 MapGeometricFields<vector, fvsPatchField, fvMeshMapper, surfaceMesh>
664 MapGeometricFields<symmTensor, fvsPatchField, fvMeshMapper, surfaceMesh>
666 MapGeometricFields<symmTensor, fvsPatchField, fvMeshMapper, surfaceMesh>
668 MapGeometricFields<tensor, fvsPatchField, fvMeshMapper, surfaceMesh>
672 MapDimensionedFields<scalar, fvMeshMapper, volMesh>(mapper);
673 MapDimensionedFields<vector, fvMeshMapper, volMesh>(mapper);
674 MapDimensionedFields<sphericalTensor, fvMeshMapper, volMesh>(mapper);
675 MapDimensionedFields<symmTensor, fvMeshMapper, volMesh>(mapper);
676 MapDimensionedFields<tensor, fvMeshMapper, volMesh>(mapper);
690 V0.setSize(nCells());
696 V0[i] = savedV0[cellMap[i]];
712 label celli = -index-2;
714 V0[celli] += savedV0[oldCelli];
722 Info<<
"Mapping old time volume V0. Merged "
723 << nMerged <<
" out of " << nCells() <<
" cells" <<
endl;
734 V00.setSize(nCells());
740 V00[i] = savedV00[cellMap[i]];
756 label celli = -index-2;
758 V00[celli] += savedV00[oldCelli];
765 Info<<
"Mapping old time volume V00. Merged "
766 << nMerged <<
" out of " << nCells() <<
" cells" <<
endl;
802 if (phiPtr_->timeIndex() != time().timeIndex())
812 scalar rDeltaT = 1.0/time().deltaTValue();
827 phibf[patchi] =
patches[patchi].patchSlice(sweptVols);
828 phibf[patchi] *= rDeltaT;
837 updateGeomNotOldVol();
841 boundary_.movePoints();
844 meshObject::movePoints<fvMesh>(*
this);
845 meshObject::movePoints<lduMesh>(*
this);
866 if (VPtr_ && (V().size() != mpm.
nOldCells()))
869 <<
"V:" << V().size()
870 <<
" not equal to the number of old cells "
874 if (V0Ptr_ && (V0Ptr_->size() != mpm.
nOldCells()))
877 <<
"V0:" << V0Ptr_->size()
878 <<
" not equal to the number of old cells "
882 if (V00Ptr_ && (V00Ptr_->size() != mpm.
nOldCells()))
885 <<
"V0:" << V00Ptr_->size()
886 <<
" not equal to the number of old cells "
897 clearGeomNotOldVol();
906 clearAddressing(
true);
908 meshObject::updateMesh<fvMesh>(*
this, mpm);
909 meshObject::updateMesh<lduMesh>(*
this, mpm);
924 ok = phiPtr_->write(valid);
932 ok = V0Ptr_->write(valid);
947 Foam::fvMesh::validComponents<Foam::sphericalTensor>()
const
int debug
Static debugging option.
The class contains the addressing required by the lduMatrix: upper, lower and losort.
void clearOut()
Clear all geometry and addressing.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual ~fvMesh()
Destructor.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
#define InfoInFunction
Report an information message using Foam::Info.
virtual bool write(const bool valid=true) const
Write mesh using IO settings from time.
void clearAddressing()
Clear topological data.
A class for managing temporary objects.
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
SolverPerformance< Type > solveSegregatedOrCoupled(const dictionary &)
Solve segregated or coupled returning the solution statistics.
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool valid) const
Write the underlying polyMesh and other data.
Selector class for finite volume differencing schemes. fvMesh is derived from fvShemes so that all fi...
Template functions to aid in the implementation of demand driven data.
bool moving() const
Is mesh moving.
virtual const fileName & dbDir() const
Override the objectRegistry dbDir for a single-region case.
virtual void mapFields(const mapPolyMesh &mpm)
Map all fields in time using given map.
label nOldFaces() const
Number of old faces.
void removeFvBoundary()
Remove boundary patches. Warning: fvPatchFields hold ref to.
bool isFile(const fileName &name, const bool checkGzip=true, const bool followLink=true)
Does the name exist as a FILE in the file system?
const fileOperation & fileHandler()
Get current file handler.
void clearOut()
Clear all geometry and addressing unnecessary for CFD.
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool operator!=(const fvMesh &rhs) const
Compares addresses.
const DimensionedField< scalar, volMesh > & V00() const
Return old-old-time cell volumes.
bool operator==(const fvMesh &rhs) const
Compares addresses.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
Generic Geometric field mapper. For "real" mapping, add template specialisations for mapping of inter...
Foam::fvMeshLduAddressing.
Registry of regIOobjects.
virtual readUpdateState readUpdate()
Update the mesh based on the mesh files saved in time.
void deleteDemandDrivenData(DataPtr &dataPtr)
static int debug
Debug switch.
virtual SolverPerformance< scalar > solve(fvMatrix< scalar > &, const dictionary &) const
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Representation of a major/minor version number.
virtual bool write(const bool valid=true) const
Write using setting from DB.
messageStream Info
Information stream (uses stdout - output is on the master only)
virtual const lduAddressing & lduAddr() const
Return ldu addressing.
virtual void updateMesh(const mapPolyMesh &mpm)
Update the mesh corresponding to given map.
Selector class for finite volume solution solution. fvMesh is derived from fvSolution so that all fie...
virtual void updateMesh(const mapPolyMesh &mpm)
Update mesh corresponding to the given map.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool valid) const
Write the objects.
streamFormat
Data format (ascii | binary)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
void addFvPatches(PtrList< polyPatch > &plist, const bool validBoundary=true)
Add boundary patches. Constructor helper.
bool movePoints()
Do what is necessary if the mesh has moved.
Mesh data needed to do the Finite Volume discretisation.
errorManip< error > abort(error &err)
const scalarField & oldCellVolumes() const
errorManipArg< error, int > exit(error &err, const int errNo=1)
void removeBoundary()
Remove boundary patches.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
readUpdateState
Enumeration defining the state of the mesh after a read update.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
const labelList & reverseCellMap() const
Reverse cell map.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Class holds all the necessary information for mapping fields associated with fvMesh.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
Traits class for primitives.
void mapClouds(const objectRegistry &db, const mapPolyMesh &mapper)
Generic Geometric field mapper.
label nOldCells() const
Number of old cells.
static void clearUpto(objectRegistry &obr)
const polyBoundaryMesh & patches
const labelList & faceMap() const
Old face map.
Cell to surface interpolation scheme. Included in fvMesh.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
compressionType
Compression treatment (UNCOMPRESSED | COMPRESSED)
const Time & time() const
Return the top-level database.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
SubField< scalar > subField
Declare type of subField.
const dimensionSet dimVolume(pow3(dimLength))
defineTypeNameAndDebug(combustionModel, 0)
Database for solution data, solver performance and other reduced data.
const labelList & cellMap() const
Old cell map.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void clearOut()
Clear all geometry and addressing.
virtual readUpdateState readUpdate()
Update the mesh based on the mesh files saved in.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
A class representing the concept of 0 (zero), which can be used to avoid manipulating objects that ar...