Go to the documentation of this file.
43 const Foam::scalar Foam::twoDPointCorrector::edgeOrthogonalityTol = 1.0 - 1
e-4;
48 void Foam::twoDPointCorrector::calcAddressing()
const
51 planeNormalPtr_ =
new vector(0, 0, 0);
52 vector& pn = *planeNormalPtr_;
65 if (isA<wedgePolyPatch>(
p))
69 const wedgePolyPatch& wp = refCast<const wedgePolyPatch>(
p);
71 pn = wp.centreNormal();
73 wedgeAxis_ = wp.axis();
74 wedgeAngle_ =
mag(
acos(wp.cosAngle()));
78 Pout<<
"Found normal from wedge patch " <<
p.index() <<
nl;
90 if (isA<emptyPolyPatch>(
p) &&
p.size())
92 pn =
p.faceAreas()[0];
96 Pout<<
"Found normal from empty patch " <<
p.index() <<
nl;
105 if (
mag(pn) < VSMALL)
108 <<
"Cannot determine normal vector from patches."
118 Pout<<
" twoDPointCorrector normal: " << pn <<
nl;
123 labelList& neIndices = *normalEdgeIndicesPtr_;
129 label nNormalEdges = 0;
133 const edge&
e = meshEdges[edgeI];
135 vector edgeVector =
e.unitVec(meshPoints);
137 if (
mag(edgeVector & pn) > edgeOrthogonalityTol)
140 neIndices[nNormalEdges++] = edgeI;
144 neIndices.setSize(nNormalEdges);
151 if (meshPoints.size() % 2 != 0)
154 <<
"the number of vertices in the geometry "
155 <<
"is odd - this should not be the case for a 2-D case. "
156 <<
"Please check the geometry."
160 if (2*nNormalEdges != meshPoints.size())
163 <<
"The number of points in the mesh is "
164 <<
"not equal to twice the number of edges normal to the plane "
165 <<
"- this may be OK only for wedge geometries.\n"
166 <<
" Please check the geometry or adjust "
167 <<
"the orthogonality tolerance.\n" <<
endl
168 <<
"Number of normal edges: " << nNormalEdges
169 <<
" number of points: " << meshPoints.size()
176 void Foam::twoDPointCorrector::clearAddressing()
const
183 void Foam::twoDPointCorrector::snapToWedge
190 scalar ADash =
mag(
A - wedgeAxis_*(wedgeAxis_ &
A));
191 vector pDash = ADash*
tan(wedgeAngle_)*planeNormal();
202 required_(mesh_.nGeometricD() == 2),
203 planeNormalPtr_(nullptr),
204 normalEdgeIndicesPtr_(nullptr),
224 const vector& pn = planeNormal();
226 if (
mag(pn.
x()) >= edgeOrthogonalityTol)
230 else if (
mag(pn.
y()) >= edgeOrthogonalityTol)
234 else if (
mag(pn.
z()) >= edgeOrthogonalityTol)
240 <<
"Plane normal not aligned with the coordinate system" <<
nl
250 if (!planeNormalPtr_)
255 return *planeNormalPtr_;
261 if (!normalEdgeIndicesPtr_)
266 return *normalEdgeIndicesPtr_;
272 if (!required_)
return;
280 const edgeList& meshEdges = mesh_.edges();
282 const labelList& neIndices = normalEdgeIndices();
283 const vector& pn = planeNormal();
285 for (
const label edgei : neIndices)
287 point& pStart =
p[meshEdges[edgei].start()];
289 point& pEnd =
p[meshEdges[edgei].end()];
292 point A = 0.5*(pStart + pEnd);
297 snapToWedge(pn,
A, pStart);
298 snapToWedge(pn,
A, pEnd);
303 pStart =
A + pn*(pn & (pStart -
A));
304 pEnd =
A + pn*(pn & (pEnd -
A));
316 if (!required_)
return;
324 const edgeList& meshEdges = mesh_.edges();
326 const labelList& neIndices = normalEdgeIndices();
327 const vector& pn = planeNormal();
329 for (
const label edgei : neIndices)
331 const edge&
e = meshEdges[edgei];
333 label startPointi =
e.start();
334 point pStart =
p[startPointi] + disp[startPointi];
336 label endPointi =
e.end();
337 point pEnd =
p[endPointi] + disp[endPointi];
340 point A = 0.5*(pStart + pEnd);
345 snapToWedge(pn,
A, pStart);
346 snapToWedge(pn,
A, pEnd);
347 disp[startPointi] = pStart -
p[startPointi];
348 disp[endPointi] = pEnd -
p[endPointi];
353 disp[startPointi] = (
A + pn*(pn & (pStart -
A))) -
p[startPointi];
354 disp[endPointi] = (
A + pn*(pn & (pEnd -
A))) -
p[endPointi];
int debug
Static debugging option.
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
virtual const pointField & points() const
Return raw points.
const Cmpt & x() const
Access to the vector x component.
dimensionedScalar tan(const dimensionedScalar &ds)
List< edge > edgeList
A List of edges.
const vector & planeNormal() const
Return plane normal.
static constexpr const zero Zero
Global zero.
void correctDisplacement(const pointField &p, vectorField &disp) const
Correct motion displacements.
~twoDPointCorrector()
Destructor.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
label nEdges() const
Number of mesh edges.
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
dimensionedScalar sign(const dimensionedScalar &ds)
const Cmpt & z() const
Access to the vector z component.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
void deleteDemandDrivenData(DataPtr &dataPtr)
Class applies a two-dimensional correction to mesh motion point field.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
errorManip< error > abort(error &err)
Vector< scalar > vector
A scalar version of the templated Vector.
const Cmpt & y() const
Access to the vector y component.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedScalar acos(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
const polyBoundaryMesh & patches
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
bool movePoints()
Correct weighting factors for moving mesh.
const labelList & normalEdgeIndices() const
Return indices of normal edges.
void correctPoints(pointField &p) const
Correct motion points.
vector point
Point is a vector.
direction normalDir() const
Return direction normal to plane.
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
void updateMesh(const mapPolyMesh &)
Update topology.