44 const Foam::scalar Foam::twoDPointCorrector::edgeOrthogonalityTol = 1.0 - 1
e-4;
49 void Foam::twoDPointCorrector::calcAddressing()
const 52 planeNormalPtr_ =
new vector(0, 0, 0);
53 vector& pn = *planeNormalPtr_;
66 if (isA<wedgePolyPatch>(
p))
70 const wedgePolyPatch& wp = refCast<const wedgePolyPatch>(
p);
72 pn = wp.centreNormal();
74 wedgeAxis_ = wp.axis();
75 wedgeAngle_ =
mag(
acos(wp.cosAngle()));
79 Pout<<
"Found normal from wedge patch " <<
p.index() <<
nl;
91 if (isA<emptyPolyPatch>(
p) &&
p.size())
93 pn =
p.faceAreas()[0];
97 Pout<<
"Found normal from empty patch " <<
p.index() <<
nl;
106 if (
mag(pn) < VSMALL)
109 <<
"Cannot determine normal vector from patches." 119 Pout<<
" twoDPointCorrector normal: " << pn <<
nl;
124 labelList& neIndices = *normalEdgeIndicesPtr_;
130 label nNormalEdges = 0;
134 const edge&
e = meshEdges[edgeI];
136 vector edgeVector =
e.unitVec(meshPoints);
138 if (
mag(edgeVector & pn) > edgeOrthogonalityTol)
141 neIndices[nNormalEdges++] = edgeI;
145 neIndices.setSize(nNormalEdges);
152 if (meshPoints.size() % 2 != 0)
155 <<
"the number of vertices in the geometry " 156 <<
"is odd - this should not be the case for a 2-D case. " 157 <<
"Please check the geometry." 161 if (2*nNormalEdges != meshPoints.size())
164 <<
"The number of points in the mesh is " 165 <<
"not equal to twice the number of edges normal to the plane " 166 <<
"- this may be OK only for wedge geometries.\n" 167 <<
" Please check the geometry or adjust " 168 <<
"the orthogonality tolerance.\n" <<
endl 169 <<
"Number of normal edges: " << nNormalEdges
170 <<
" number of points: " << meshPoints.size()
177 void Foam::twoDPointCorrector::clearAddressing()
const 184 void Foam::twoDPointCorrector::snapToWedge
191 scalar ADash =
mag(
A - wedgeAxis_*(wedgeAxis_ &
A));
192 vector pDash = ADash*
tan(wedgeAngle_)*planeNormal();
200 Foam::twoDPointCorrector::twoDPointCorrector(
const polyMesh&
mesh)
203 required_(mesh_.nGeometricD() == 2),
204 planeNormalPtr_(nullptr),
205 normalEdgeIndicesPtr_(nullptr),
225 const vector& pn = planeNormal();
227 if (
mag(pn.x()) >= edgeOrthogonalityTol)
231 else if (
mag(pn.y()) >= edgeOrthogonalityTol)
235 else if (
mag(pn.z()) >= edgeOrthogonalityTol)
241 <<
"Plane normal not aligned with the coordinate system" <<
nl 251 if (!planeNormalPtr_)
256 return *planeNormalPtr_;
262 if (!normalEdgeIndicesPtr_)
267 return *normalEdgeIndicesPtr_;
273 if (!required_)
return;
281 const edgeList& meshEdges = mesh_.edges();
283 const labelList& neIndices = normalEdgeIndices();
284 const vector& pn = planeNormal();
286 for (
const label edgei : neIndices)
288 point& pStart =
p[meshEdges[edgei].start()];
290 point& pEnd =
p[meshEdges[edgei].end()];
293 point A = 0.5*(pStart + pEnd);
298 snapToWedge(pn,
A, pStart);
299 snapToWedge(pn,
A, pEnd);
304 pStart =
A + pn*(pn & (pStart -
A));
305 pEnd =
A + pn*(pn & (pEnd -
A));
317 if (!required_)
return;
325 const edgeList& meshEdges = mesh_.edges();
327 const labelList& neIndices = normalEdgeIndices();
328 const vector& pn = planeNormal();
330 for (
const label edgei : neIndices)
332 const edge&
e = meshEdges[edgei];
334 label startPointi =
e.start();
335 point pStart =
p[startPointi] + disp[startPointi];
337 label endPointi =
e.end();
338 point pEnd =
p[endPointi] + disp[endPointi];
341 point A = 0.5*(pStart + pEnd);
346 snapToWedge(pn,
A, pStart);
347 snapToWedge(pn,
A, pEnd);
348 disp[startPointi] = pStart -
p[startPointi];
349 disp[endPointi] = pEnd -
p[endPointi];
354 disp[startPointi] = (
A + pn*(pn & (pStart -
A))) -
p[startPointi];
355 disp[endPointi] = (
A + pn*(pn & (pEnd -
A))) -
p[endPointi];
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
dimensionedScalar sign(const dimensionedScalar &ds)
void updateMesh(const mapPolyMesh &)
Update topology.
dimensionedScalar acos(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr char nl
The newline '\n' character (0x0a)
const labelList & normalEdgeIndices() const
Return indices of normal edges.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
virtual const pointField & points() const
Return raw points.
#define forAll(list, i)
Loop across all elements in list.
void correctPoints(pointField &p) const
Correct motion points.
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
void correctDisplacement(const pointField &p, vectorField &disp) const
Correct motion displacements.
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
List< edge > edgeList
A List of edges.
Class applies a two-dimensional correction to mesh motion point field.
errorManip< error > abort(error &err)
label nEdges() const
Number of mesh edges.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
const vector & planeNormal() const
Return plane normal.
defineTypeNameAndDebug(combustionModel, 0)
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
Template functions to aid in the implementation of demand driven data.
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.
bool movePoints()
Correct weighting factors for moving mesh.
PtrList< volScalarField > & Y
direction normalDir() const
Return direction normal to plane.
const polyBoundaryMesh & patches
Field< vector > vectorField
Specialisation of Field<T> for vector.
Mesh consisting of general polyhedral cells.
~twoDPointCorrector()
Destructor.
List< label > labelList
A List of labels.
dimensionedScalar tan(const dimensionedScalar &ds)
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
void deleteDemandDrivenData(DataPtr &dataPtr)
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
static constexpr const zero Zero
Global zero (0)