35 void Foam::plane::makeUnitNormal
37 const char *
const caller,
41 const scalar magNormal(
Foam::mag(normal_));
43 if (magNormal < VSMALL)
46 <<
"Plane normal has zero length.\nCalled from " << caller
57 void Foam::plane::calcFromCoeffs
63 const char *
const caller
68 origin_ =
vector((-d/a), 0, 0);
70 else if (
mag(
b) > VSMALL)
74 else if (
mag(
c) > VSMALL)
81 <<
"At least one plane coefficient must have a value" 86 makeUnitNormal(caller);
90 void Foam::plane::calcFromEmbeddedPoints
95 const char *
const caller
98 origin_ = (point1 + point2 + point3)/3;
99 const vector line12 = point1 - point2;
100 const vector line23 = point2 - point3;
105 ||
mag(line23) < VSMALL
106 ||
mag(point3-point1) < VSMALL
110 <<
"Bad points:" << point1 <<
' ' << point2 <<
' ' << point3
114 normal_ = line12 ^ line23;
116 makeUnitNormal(caller);
124 normal_(normalVector),
133 const point& originPoint,
134 const vector& normalVector,
135 const bool doNormalise
138 normal_(normalVector),
183 dict.readIfPresent(
"planeType", planeType);
185 if (planeType.empty())
187 const dictionary& coeffs =
dict.optionalSubDict(
"pointAndNormalDict");
189 origin_ = coeffs.
get<
point>(
"point");
190 normal_ = coeffs.
get<
point>(
"normal");
192 makeUnitNormal(
"point/normal");
194 else if (planeType ==
"pointAndNormal")
196 const dictionary& coeffs =
dict.subDict(
"pointAndNormalDict");
198 origin_ = coeffs.getCompat<
point>(
"point", {{
"basePoint", 1612}});
199 normal_ = coeffs.getCompat<
point>(
"normal", {{
"normalVector", 1612}});
201 makeUnitNormal(
"point/normal");
203 else if (planeType ==
"planeEquation")
205 const dictionary& subDict =
dict.subDict(
"planeEquationDict");
209 subDict.get<scalar>(
"a"),
210 subDict.get<scalar>(
"b"),
211 subDict.get<scalar>(
"c"),
212 subDict.get<scalar>(
"d"),
216 else if (planeType ==
"embeddedPoints")
218 const dictionary& subDict =
dict.subDict(
"embeddedPointsDict");
220 calcFromEmbeddedPoints
222 subDict.get<
point>(
"point1"),
223 subDict.get<
point>(
"point2"),
224 subDict.get<
point>(
"point3"),
231 <<
"Invalid plane type: " << planeType <<
nl 232 <<
"Valid options: (planeEquation embeddedPoints pointAndNormal)" 253 const scalar magX =
mag(normal_.x());
254 const scalar magY =
mag(normal_.y());
255 const scalar magZ =
mag(normal_.z());
262 coeffs[1] = normal_.y()/normal_.x();
263 coeffs[2] = normal_.z()/normal_.x();
267 coeffs[0] = normal_.x()/normal_.z();
268 coeffs[1] = normal_.y()/normal_.z();
276 coeffs[0] = normal_.x()/normal_.y();
278 coeffs[2] = normal_.z()/normal_.y();
282 coeffs[0] = normal_.x()/normal_.z();
283 coeffs[1] = normal_.y()/normal_.z();
290 - coeffs[0] * origin_.x()
291 - coeffs[1] * origin_.y()
292 - coeffs[2] * origin_.z()
305 const scalar denom =
stabilise((dir & normal_), VSMALL);
307 return ((origin_ - pnt0) & normal_)/denom;
318 const vector& n1 = this->normal();
321 const point& p1 = this->origin();
324 const scalar n1p1 = n1 & p1;
325 const scalar n2p2 = n2 & p2;
327 const vector dir = n1 ^ n2;
331 const scalar magX =
mag(dir.x());
332 const scalar magY =
mag(dir.y());
333 const scalar magZ =
mag(dir.z());
372 pt[i1] = (n2[i2]*n1p1 - n1[i2]*n2p2) / (n1[i1]*n2[i2] - n2[i1]*n1[i2]);
373 pt[i2] = (n2[i1]*n1p1 - n1[i1]*n2p2) / (n1[i2]*n2[i1] - n1[i1]*n2[i2]);
391 coeffs1[0],coeffs1[1],coeffs1[2],
392 coeffs2[0],coeffs2[1],coeffs2[2],
393 coeffs3[0],coeffs3[1],coeffs3[2]
396 vector b(coeffs1[3],coeffs2[3],coeffs3[3]);
398 return (
inv(a) & (-
b));
408 point p = origin_ + point::uniform(dist);
414 p[0] = -1.0*(coeff[1]*
p[1] + coeff[2]*
p[2] + coeff[3])/coeff[0];
418 p[1] = -1.0*(coeff[0]*
p[0] + coeff[2]*
p[2] + coeff[3])/coeff[1];
423 p[2] = -1.0*(coeff[0]*
p[0] + coeff[1]*
p[1] + coeff[3])/coeff[2];
432 const vector mirroredPtDir =
p - nearestPoint(
p);
434 if ((normal() & mirroredPtDir) > 0)
A reference point and direction.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A 1D vector of objects of type <T> with a fixed length <N>.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
constexpr char nl
The newline '\n' character (0x0a)
dimensionedScalar stabilise(const dimensionedScalar &x, const dimensionedScalar &y)
scalar distance(const vector &p1, const vector &p2)
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Geometric class that creates a 3D plane and can return the intersection point between a line and the ...
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T. FatalIOError if not found, or if the number of tokens is incorrect.
const point & origin() const noexcept
The plane base point.
scalar normalIntersect(const point &pnt0, const vector &dir) const
Return cut coefficient for plane and line defined by origin and direction.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
A class for handling words, derived from Foam::string.
virtual Ostream & endBlock()
Write end block group.
void writeDict(Ostream &os) const
Write to dictionary.
errorManip< error > abort(error &err)
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
plane()
Construct zero-initialised.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
OBJstream os(runTime.globalPath()/outputName)
point mirror(const point &p) const
Mirror the supplied point in the plane. Return the mirrored point.
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
vector point
Point is a vector.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
const dimensionedScalar c
Speed of light in a vacuum.
virtual Ostream & beginBlock(const keyType &kw)
Write begin block group with the given name.
ray planeIntersect(const plane &plane2) const
Return the cutting line between this plane and another.
point somePointInPlane(const scalar dist=1e-3) const
Return a point somewhere on the plane, a distance from the base.
FixedList< scalar, 4 > planeCoeffs() const
Return coefficients for the plane equation: ax + by + cz + d = 0.
point planePlaneIntersect(const plane &plane2, const plane &plane3) const
Return the cutting point between this plane and two other planes.
Tensor of scalars, i.e. Tensor<scalar>.
const vector & normal() const noexcept
The plane unit normal.
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...
static constexpr const zero Zero
Global zero (0)