Go to the documentation of this file.
60 { outputFormatType::PLAIN,
"plain" },
61 { outputFormatType::DICTIONARY,
"dictionary" },
71 { scalingType::LENGTH,
"length" },
72 { scalingType::FORCE,
"force" },
73 { scalingType::MOMENT,
"moment" },
85 static inline Ostream& putPlain(Ostream& os,
const vector& v)
87 return os << v.x() <<
' ' << v.y() <<
' ' << v.z();
95 static void writeList(Ostream& os,
const string& header,
const UList<T>& list)
97 const label len = list.size();
100 os << header.c_str() <<
nl;
106 for (label i=0; i < len; ++i)
133 inputName_(
"positions.in"),
134 outputName_(
"forces.out"),
135 logName_(
"movement.log"),
161 inputName_(
"positions.in"),
162 outputName_(
"forces.out"),
163 logName_(
"movement.log"),
177 return (
timeIndex >= lastTrigger_ + calcFrequency_);
193 auto&
points0 = tpoints0.ref();
198 originalIds_.clear();
199 controllers_.clear();
200 patchControls_.clear();
208 if (originalIds_.size() !=
points0.size())
211 <<
"Incorrect number of pointLabels. Had "
212 << originalIds_.size() <<
" for " <<
points0.size() <<
" points"
221 for (
const label
id : originalIds_)
223 if (!pointIdMap.insert(
id, pointi))
231 if (!duplicates.
empty())
234 <<
"Found duplicate point ids "
248 (*iter)->remapPointLabels(
points0.size(), pointIdMap);
256 controllers_.clear();
265 coupler_.readDict(commDict);
267 calcFrequency_ = commDict.
getOrDefault<label>(
"calcFrequency", 1);
271 commDict.
readEntry(
"inputName", inputName_);
272 commDict.
readEntry(
"outputName", outputName_);
276 outputFormat_ = formatNames.get(
"outputFormat", commDict);
283 if ((scaleDict = commDict.
findDict(
"scaleInput")) !=
nullptr)
285 for (
int i=0; i < scaleInput_.size(); ++i)
292 && scaleInput_[i] > 0
295 Info<<
"Using input " << key <<
" multiplier: "
296 << scaleInput_[i] <<
nl;
301 if ((scaleDict = commDict.
findDict(
"scaleOutput")) !=
nullptr)
303 for (
int i=0; i < scaleOutput_.size(); ++i)
310 && scaleOutput_[i] > 0
313 Info<<
"Using output " << key <<
" multiplier: "
314 << scaleOutput_[i] <<
nl;
326 quaternion::eulerOrder::ZXZ
331 state0_.
scalePoints(scaleInput_[scalingType::LENGTH]);
342 return hasPatchControl(pp.
index());
351 return hasInterpolator(fpatch.
index());
360 const auto ctrlIter = patchControls_.cfind(pp.
index());
362 if (!ctrlIter.good())
365 <<
"No controllers for patch " << pp.
name()
369 const patchControl& ctrl = *ctrlIter;
371 for (
const word& ctrlName : ctrl.names_)
373 const auto iter = controllers_.cfind(ctrlName);
378 <<
"No controller: " << ctrlName <<
nl
379 <<
" For patch " << pp.
name()
394 const pointField& lumpedCentres0 = state0().points();
396 const label patchIndex = pp.
index();
401 patchControl& ctrl = patchControls_(patchIndex);
402 ctrl.names_ = ctrlNames;
404 labelList& faceToPoint = ctrl.faceToPoint_;
405 faceToPoint.
resize(pp.size(), -1);
407 checkPatchControl(pp);
414 for (
const word& ctrlName : ctrl.names_)
416 const auto iter = controllers_.
cfind(ctrlName);
421 <<
"No controller: " << ctrlName <<
nl
430 if (ctrl.names_.size() && subsetPointIds.
empty())
433 <<
"Controllers specified, but without any points" <<
nl
442 subsetPointIds.
size()
458 const scalar searchDistSqr(
sqr(GREAT));
460 const label patchStart = pp.
start();
463 const point fc(faces[patchStart + patchFacei].centre(
points0));
466 faceToPoint[patchFacei] =
467 treePoints.pointLabel
469 ppTree.findNearest(fc, searchDistSqr).index()
475 Pout<<
"Added face mapping for patch: " << patchIndex <<
endl;
487 const pointField& lumpedCentres0 = state0().points();
489 const label patchIndex = fpatch.
index();
491 patchControl& ctrl = patchControls_(patchIndex);
502 for (
const word& ctrlName : ctrl.names_)
504 const auto iter = controllers_.
cfind(ctrlName);
509 <<
"No controller: " << ctrlName <<
nl
525 adjacency(prev).insert(curr);
526 adjacency(curr).insert(prev);
528 else if (!adjacency.found(curr))
530 adjacency(curr).clear();
535 if (ctrl.names_.empty())
538 const label len = state0().size();
540 for (label i=0; i < len; ++i)
542 const label curr = i;
546 const label prev = i-1;
547 adjacency(prev).insert(curr);
548 adjacency(curr).insert(prev);
550 else if (!adjacency.found(curr))
552 adjacency(curr).clear();
557 if (ctrl.names_.size() && subsetPointIds.
empty())
560 <<
"Controllers specified, but without any points" <<
nl
580 const label nearest = iter.key();
582 labelList neighbours(iter.val().sortedToc());
584 const label len = neighbours.size();
589 const point& nearPt = lumpedCentres0[nearest];
591 for (label j=1; j < len; ++j)
593 for (label i=0; i < j; ++i)
595 labelPair neiPair(neighbours[i], neighbours[j]);
600 lumpedCentres0[neiPair.first()],
601 lumpedCentres0[neiPair.
second()]
605 if (tri.pointToBarycentric(tri.a(), bary) > SMALL)
616 acuteAngles.
append(-(ab & ac));
621 if (pairs.size() > 1)
628 adjacencyPairs.insert(nearest, pairs);
634 Info<<
"Adjacency table for patch: " << fpatch.
name() <<
nl;
636 for (
const label own : adjacency.sortedToc())
639 for (
const label nei : adjacency[own].sortedToc())
644 if (originalIds_.size())
646 Info<<
" # " << originalIds_[own] <<
" =>";
647 for (
const label nei : adjacency[own].sortedToc())
649 Info<<
' ' << originalIds_[nei];
662 subsetPointIds.
size()
680 const scalar searchDistSqr(
sqr(GREAT));
684 interpList.
resize(meshPoints.size());
689 forAll(meshPoints, pointi)
694 const label nearest =
695 treePoints.pointLabel
697 ppTree.findNearest(ptOnMesh, searchDistSqr).index()
700 interpList[pointi].nearest(nearest);
717 const point& nearPt = lumpedCentres0[nearest];
721 const labelPairList& adjacentPairs = adjacencyPairs[nearest];
723 unsortedNeighbours = adjacency[nearest].toc();
724 unsortedNeiWeightDist.
resize(unsortedNeighbours.size());
726 forAll(unsortedNeighbours, nbri)
728 unsortedNeiWeightDist[nbri] =
729 magSqr(ptOnMesh - lumpedCentres0[unsortedNeighbours[nbri]]);
738 for (
const label nbri : distOrder)
740 const label nextPointi = unsortedNeighbours[nbri];
742 const point& nextPt = lumpedCentres0[nextPointi];
746 const scalar weight =
757 unsortedNeiWeightDist[nbri] = weight;
760 distOrder[ngood] = nbri;
767 if (distOrder.size() < 1)
775 bool useFirst =
true;
777 if (neighbours.size() > 1 && adjacentPairs.size())
782 neiPointid.
set(neighbours);
784 for (
const labelPair& ends : adjacentPairs)
786 label nei1 = ends.first();
787 label nei2 = ends.second();
789 if (!neiPointid.
test(nei1) || !neiPointid.
test(nei2))
794 else if (neighbours.find(nei2) < neighbours.find(nei1))
802 triFace triF(nearest, nei1, nei2);
812 interpList[pointi].
set(triF, bary);
823 interpList[pointi].next(neighbours.first(), neiWeight.first());
834 const label nLumpedPoints = state0().size();
838 if (patchControls_.empty())
841 <<
"Attempted to calculate areas without setMapping()"
849 if (isA<fvMesh>(pmesh))
851 const fvMesh&
mesh = dynamicCast<const fvMesh>(pmesh);
854 const surfaceVectorField::Boundary& patchSf =
859 const label patchIndex = iter.key();
860 const patchControl& ctrl = iter.val();
862 const labelList& faceToPoint = ctrl.faceToPoint_;
869 const label pointIndex = faceToPoint[patchFacei];
878 zoneAreas[pointIndex] +=
mag(patchSf[patchIndex][patchFacei]);
897 const label nLumpedPoints = state0().size();
899 forces.
resize(nLumpedPoints);
900 moments.
resize(nLumpedPoints);
902 if (patchControls_.empty())
905 <<
"Attempted to calculate forces without setMapping()"
922 const word pName(forcesDict_.getOrDefault<
word>(
"p",
"p"));
923 scalar pRef = forcesDict_.getOrDefault<scalar>(
"pRef", 0);
924 scalar rhoRef = forcesDict_.getOrDefault<scalar>(
"rhoRef", 1);
933 if (isA<fvMesh>(pmesh) && pPtr)
935 const fvMesh&
mesh = dynamicCast<const fvMesh>(pmesh);
948 rhoRef = (
p.dimensions() ==
dimPressure ? 1.0 : rhoRef);
955 const label patchIndex = iter.key();
956 const patchControl& ctrl = iter.val();
958 const labelList& faceToPoint = ctrl.faceToPoint_;
960 if (!forceOnPatches.found(patchIndex))
970 * patchSf[patchIndex] * (patchPress[patchIndex] - pRef)
975 const vectorField& forceOnPatch = *forceOnPatches[patchIndex];
982 const label pointIndex = faceToPoint[patchFacei];
991 forces[pointIndex] += forceOnPatch[patchFacei];
998 patchCf[patchIndex][patchFacei]
999 - lumpedCentres[pointIndex]
1003 moments[pointIndex] += lever ^ forceOnPatch[patchFacei];
1029 return pointsDisplacement(state(), fpatch,
points0);
1041 const label patchIndex = fpatch.
index();
1044 const pointField& lumpedCentres0 = state0().points();
1057 auto& disp = tdisp.ref();
1061 = patchControls_[patchIndex].interp_;
1063 forAll(meshPoints, pointi)
1073 disp[pointi] = (rotTensor & (
p0 - origin0)) + origin -
p0;
1088 const label patchIndex = fpatch.
index();
1091 const pointField& lumpedCentres0 = state0().points();
1104 auto& disp = tdisp.ref();
1108 patchControls_[patchIndex].interp_;
1110 forAll(meshPoints, pointi)
1120 disp[pointi] = (rotTensor & (
p0 - origin0)) + origin;
1138 const bool status = state_.
readData
1141 coupler().resolveFile(inputName_),
1142 state0().rotationOrder(),
1146 scalePoints(state_);
1148 state_.relax(relax_, prev);
1160 const Time* timeinfo
1163 const bool writeMoments = (moments.
size() == forces.
size());
1165 if (fmt == outputFormatType::PLAIN)
1167 os <<
"########" <<
nl;
1170 const Time& t = *timeinfo;
1175 os <<
"# size=" << this->size() <<
nl
1176 <<
"# columns (points) (forces)";
1185 bool report =
false;
1186 scalar scaleLength = scaleOutput_[scalingType::LENGTH];
1187 scalar scaleForce = scaleOutput_[scalingType::FORCE];
1188 scalar scaleMoment = scaleOutput_[scalingType::MOMENT];
1190 if (scaleLength > 0)
1210 if (scaleMoment > 0)
1222 os <<
"# scaling points=" << scaleLength
1223 <<
" forces=" << scaleForce;
1227 os <<
" moments=" << scaleMoment;
1233 os <<
"########" <<
nl;
1237 const point& pt = state0().points()[i];
1239 putPlain(os, scaleLength * pt) <<
' ';
1241 if (i < forces.
size())
1243 const vector val(scaleForce * forces[i]);
1254 if (i < moments.
size())
1256 const vector val(scaleMoment * moments[i]);
1274 os <<
"////////" <<
nl;
1277 const Time& t = *timeinfo;
1301 const Time* timeinfo
1311 const fileName output(coupler().resolveFile(outputName_));
1314 writeData(os, forces, moments, outputFormat_, timeinfo);
1319 const fileName output(coupler().resolveFile(logName_));
1330 writeData(os, forces, moments, outputFormatType::PLAIN, timeinfo);
int debug
Static debugging option.
const T & second() const noexcept
Return second element, which is also the last element.
dictionary * findDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX)
Find and return a sub-dictionary pointer if present.
scalar timeOutputValue() const
Return current time value.
label size() const noexcept
The number of elements in table.
const dimensionSet dimPressure
virtual label index() const =0
Return the index of this patch in the pointBoundaryMesh.
void writeDict(Ostream &os) const
Write axis, locations, division as a dictionary.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
A class for handling file names.
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
bool forcesAndMoments(const polyMesh &pmesh, List< vector > &forces, List< vector > &moments) const
The forces and moments acting on each pressure-zone.
static const Enum< scalingType > scalingNames
Names for the scaling types.
T interpolate(const UList< T > &input) const
Linear interpolated value between nearest and next locations.
void readDict(const dictionary &dict)
Update settings from dictionary.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
int debugSwitch(const char *name, const int deflt=0)
Lookup debug switch or add default value.
Standard boundBox with extra functionality for use in octree.
List< scalar > areas(const polyMesh &pmesh) const
The areas for each pressure-zone.
static const versionNumber currentVersion
The current version number (2.0)
Templated 2D Barycentric derived from VectorSpace. Has 3 components, one of which is redundant.
void set(const bitSet &bitset)
Set specified bits from another bitset.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
const pointField & points() const
The points corresponding to mass centres.
Ostream & endl(Ostream &os)
Add newline and flush stream.
scalingType
Output format types.
Basic pointPatch represents a set of points from the mesh.
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
EnumType get(const word &enumName) const
The enumeration corresponding to the given name.
bool test(const label pos) const
Test value at specified position, never auto-vivify entries.
void Swap(DynamicList< T, SizeMin1 > &a, DynamicList< T, SizeMin2 > &b)
triPointRef tri(const UList< point > &points) const
Return the triangle.
Mesh consisting of general polyhedral cells.
Holds (reference to) pointField. Encapsulation of data needed for octree searches....
#define forAll(list, i)
Loop across all elements in list.
bool hasInterpolator(const pointPatch &fpatch) const
Check if patch control exists for specified patch.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
bool readData(Istream &is, const quaternion::eulerOrder rotOrder=quaternion::eulerOrder::ZXZ, const bool degrees=false)
Read input as dictionary content.
bool outside() const
True if any coordinates are negative.
bool couplingPending(const label timeIndex) const
Check if coupling is pending (according to the calcFrequency)
A triangle primitive used to calculate face normals and swept volumes.
The state of lumped points corresponds to positions and rotations.
void scalePoints(const scalar scaleFactor)
Scale points by given factor.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
static const Enum< inputFormatType > formatNames
Names for the input format types.
messageStream Info
Information stream (uses stdout - output is on the master only)
static int debug
Debug switch.
A patch is a list of labels that address the faces in the global face list.
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
void clear()
Clear the addressed list, i.e. set the size to zero.
static MinMax< T > zero_one()
A 0-1 range corresponding to the pTraits zero, one.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Non-pointer based hierarchical recursive searching.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
const polyMesh & mesh() const
Return the mesh reference.
void writeList(vtk::formatter &fmt, const UList< uint8_t > &values)
Write a list of uint8_t values.
scalar mag() const
Return scalar magnitude.
Point unitVec() const
Return the unit vector (start-to-end)
void resize(const label newSize)
Adjust allocated size of list.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
bool writeData(Ostream &os, const UList< vector > &forces, const UList< vector > &moments, const outputFormatType fmt=outputFormatType::PLAIN, const Time *timeinfo=nullptr) const
Write points, forces, moments. Only call from the master process.
static void writeData(Ostream &os, const Type &val)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
static const word canonicalName
The canonical name ("lumpedPointMovement") for the dictionary.
static const Enum< outputFormatType > formatNames
Names for the output format types.
bool readState()
Read state from file, applying relaxation as requested.
Mesh data needed to do the Finite Volume discretisation.
virtual const word & name() const =0
Return name.
Vector< scalar > vector
A scalar version of the templated Vector.
List< Key > sortedToc() const
The table of contents (the keys) in sorted order.
label start() const
Return start label of this patch in the polyMesh face list.
Inter-processor communication reduction functions.
dictionary subOrEmptyDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX, const bool mandatory=false) const
errorManipArg< error, int > exit(error &err, const int errNo=1)
void resize(const label nElem)
Alter addressable list size.
Output to file stream, using an OSstream.
virtual label size() const =0
Return size.
tmp< pointField > pointsPosition(const lumpedPointState &state, const pointPatch &fpatch, const pointField &points0) const
The points absolute position according to specified state.
static bool master(const label communicator=0)
Am I the master process.
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
pointField points0(pointIOField(IOobject("points", mesh.time().constant(), polyMesh::meshSubDir, mesh, IOobject::MUST_READ, IOobject::NO_WRITE, false)))
static const Enum< eulerOrder > eulerOrderNames
The names for Euler-angle rotation order.
virtual const faceList & faces() const
Return raw faces.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const_iterator cfind(const Key &key) const
Find and return an const_iterator set at the hashed entry.
const Type * findObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
forAllConstIters(mixture.phases(), phase)
A triangular face using a FixedList of labels corresponding to mesh vertices.
Input/output from file streams.
FlatOutput< Container > flatOutput(const Container &obj, label len=0)
Global flatOutput function.
A HashTable of pointers to objects of type <T>.
outputFormatType
Output format types.
bool hasPatchControl(const label patchIndex) const
Check if patch control exists for specified patch.
virtual const labelList & meshPoints() const =0
Return mesh points.
void setPatchControl(const polyPatch &pp, const wordList &ctrlNames, const pointField &points0)
Define pressure-zones mapping for faces in the specified patches.
tmp< pointField > pointsDisplacement(const pointPatch &fpatch, const pointField &points0) const
Displace points according to the current state.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
lumpedPointMovement()
Default construct.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
void clear()
Clear the list, i.e. set size to zero.
const tensorField & rotations() const
The local-to-global transformation for each point.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
const polyBoundaryMesh & patches
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
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 size(const label n) noexcept
Override size to be inconsistent with allocated storage.
static const Vector< scalar > zero
A List with indirect addressing.
A simple linear interpolator between two locations, which are referenced by index.
void inplaceReorder(const labelUList &oldToNew, ListType &input, const bool prune=false)
Inplace reorder the elements of a list.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
bool empty() const noexcept
True if the hash table is empty.
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
barycentric2D pointToBarycentric(const point &pt) const
Calculate the barycentric coordinates from the given point.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
void checkPatchControl(const polyPatch &pp) const
Check if patch control exists for specified patch.
const volScalarField & p0
void couplingCompleted(const label timeIndex) const
Register that coupling is completed at this calcFrequency.
label timeIndex() const
Return current time index.
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Begin list [isseparator].
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
void setInterpolator(const pointPatch &fpatch, const pointField &points0)
Check if patch control exists for specified patch.
const word & name() const
The patch name.
Point vec() const
Return start-to-end vector.
#define WarningInFunction
Report a warning using Foam::Warning.
label index() const
The index of this patch in the boundaryMesh.
labelList pointLabels(nPoints, -1)
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
bool set(const Key &key)
Same as insert (no value to overwrite)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A HashTable of pointers to objects of type <T> with a label key.
const surfaceVectorField & Sf() const
Return cell face area vectors.