48void Foam::solidBodyFvGeometryScheme::setMeshMotionData()
50 if (!cacheInitialised_ || !cacheMotion_)
54 changedFaceIDs_.
clear();
55 changedPatchIDs_.
clear();
56 changedCellIDs_.
clear();
61 if (oldPoints.size() != currPoints.size())
64 <<
"Old and current points sizes must be the same. "
65 <<
"Old points:" << oldPoints.size()
66 <<
" Current points:" << currPoints.size()
70 bitSet changedPoints(oldPoints.size());
73 forAll(changedPoints, pointi)
75 changedPoints.set(pointi, oldPoints[pointi] != currPoints[pointi]);
79 <<
"SBM --- Changed points:"
97 for (
const label pointi : changedPoints)
99 for (
const auto facei : pointFaces[pointi])
111 changedCellIDs_ =
cellIDs.toc();
114 <<
"SBM --- Changed cells:"
121 const auto changedFaceFlag = faceIDs.values();
123 DynamicList<label> changedFaceIDs(faceIDs.count());
124 DynamicList<label> changedPatchIDs(faceIDs.count());
127 if (changedFaceFlag[facei])
129 changedFaceIDs.append(facei);
130 changedPatchIDs.append(-1);
135 for (label patchi = 0; patchi < pbm.size(); ++patchi)
137 const polyPatch& pp = pbm[patchi];
139 for (
const label meshFacei : pp.range())
141 if (changedFaceFlag[meshFacei])
143 changedFaceIDs.
append(meshFacei);
144 changedPatchIDs.append(patchi);
149 changedFaceIDs_.
transfer(changedFaceIDs);
150 changedPatchIDs_.
transfer(changedPatchIDs);
153 cacheInitialised_ =
true;
166 partialUpdate_(
dict.getOrDefault<
bool>(
"partialUpdate", true)),
167 cacheMotion_(
dict.getOrDefault<
bool>(
"cacheMotion", true)),
168 cacheInitialised_(false),
174 <<
"partialUpdate:" << partialUpdate_
175 <<
" cacheMotion:" << cacheMotion_
188 mesh_.hasCellCentres()
189 && mesh_.hasFaceCentres()
190 && mesh_.hasCellVolumes()
191 && mesh_.hasFaceAreas();
196 <<
"Creating initial geometry using primitiveMesh::updateGeom"
209 const pointField& oldPoints = mesh_.oldPoints();
210 const pointField& currPoints = mesh_.points();
212 if (oldPoints.
size() != currPoints.
size())
215 <<
"Old and current points sizes must be the same. "
216 <<
"Old points:" << oldPoints.
size()
217 <<
" Current points:" << currPoints.
size()
222 const faceList& faces = mesh_.faces();
224 auto tmeshPhi(
const_cast<fvMesh&
>(mesh_).setPhi());
227 const scalar rdt = 1.0/mesh_.time().deltaTValue();
230 auto& meshPhi = tmeshPhi.ref();
231 auto& meshPhii = meshPhi.primitiveFieldRef();
232 auto& meshPhiBf = meshPhi.boundaryFieldRef();
238 forAll(changedFaceIDs_, i)
240 const face&
f = faces[changedFaceIDs_[i]];
242 if (changedPatchIDs_[i] == -1)
244 const label facei = changedFaceIDs_[i];
245 meshPhii[facei] =
f.sweptVol(oldPoints, currPoints)*rdt;
249 const label patchi = changedPatchIDs_[i];
252 if (isA<emptyPolyPatch>(pp))
257 const label patchFacei = changedFaceIDs_[i] - pp.
start();
259 meshPhiBf[patchi][patchFacei] =
260 f.sweptVol(oldPoints, currPoints)*rdt;
265 if (partialUpdate_ && haveGeometry)
300 std::move(faceCentres),
301 std::move(faceAreas),
302 std::move(cellCentres),
303 std::move(cellVolumes)
308 for (
const auto&
p : mesh_.boundaryMesh())
311 <<
" sum(Sf)=" <<
sum(
p.faceAreas())
312 <<
" sum(mag(Sf))=" <<
sum(
mag(
p.faceAreas()))
320 <<
"Performing complete geometry clear and update" <<
endl;
342 cacheInitialised_ =
false;
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
const word & name() const noexcept
Return the object name.
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
void transfer(List< T > &list)
void clear()
Clear the list, i.e. set size to zero.
void append(T *ptr)
Append an element to the end of the list.
void size(const label n)
Older name for setAddressableSize.
Default geometry calculation scheme. Slight stabilisation for bad meshes.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
A face is a list of labels corresponding to mesh vertices.
Abstract base class for geometry calculation schemes.
const fvMesh & mesh_
Hold reference to mesh.
Mesh data needed to do the Finite Volume discretisation.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
void updateMesh()
Update for new mesh topology.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
virtual const labelList & faceOwner() const
Return face owner.
virtual const pointField & oldPoints() const
Return old points (mesh motion)
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const labelList & faceNeighbour() const
Return face neighbour.
virtual const pointField & points() const
Return raw points.
A patch is a list of labels that address the faces in the global face list.
label start() const
Return start label of this patch in the polyMesh face list.
void clearGeom()
Clear geometry.
label nInternalFaces() const noexcept
Number of internal faces.
label nCells() const noexcept
Number of mesh cells.
label nFaces() const noexcept
Number of mesh faces.
const labelListList & pointFaces() const
virtual void updateGeom()
Update all geometric data.
void resetGeometry(pointField &&faceCentres, pointField &&faceAreas, pointField &&cellCentres, scalarField &&cellVolumes)
Reset the local geometry.
Geometry calculation scheme that performs geometry updates only in regions where the mesh has changed...
virtual void movePoints()
Do what is necessary if the mesh has moved.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
#define DebugInFunction
Report an information message using Foam::Info.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
vectorField pointField
pointField is a vectorField.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManip< error > abort(error &err)
static constexpr const zero Zero
Global zero (0)
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
#define forAll(list, i)
Loop across all elements in list.