47 void Foam::blockEdges::arcEdge::calcFromMidPoint
58 const scalar asqr = a & a;
59 const scalar bsqr =
b &
b;
60 const scalar adotb = a &
b;
62 const scalar denom = asqr*bsqr - adotb*adotb;
64 if (
mag(denom) < ROOTVSMALL)
71 const scalar fact = 0.5*(bsqr - adotb)/denom;
73 const point centre = p1 + 0.5*a + fact*((a ^
b) ^ a);
76 const vector r1(p1 - centre);
77 const vector r2(p2 - centre);
78 const vector r3(p3 - centre);
80 const scalar mag1(
mag(r1));
81 const scalar mag3(
mag(r3));
90 angle_ =
acos((r1 & r3)/(mag1*mag3));
93 if (((r1 ^ r2) & (r1 ^ r3)) < 0.0)
100 if (
mag(arcAxis)/(mag1*mag3) < 0.001)
111 cs_ = coordSystem::cylindrical(centre, arcAxis, r1);
115 void Foam::blockEdges::arcEdge::calcFromCentre
125 const vector r1(p1 - centre);
126 const vector r3(p3 - centre);
128 const scalar mag1(
mag(r1));
129 const scalar mag3(
mag(r3));
131 const vector chord(p3 - p1);
133 const vector arcAxis(r1 ^ r3);
136 radius_ = 0.5*(mag1 + mag3);
139 angle_ =
acos((r1 & r3)/(mag1*mag3));
143 bool needsAdjust =
false;
147 needsAdjust = !
equal(mag1, mag3);
149 if (!
equal(rMultiplier, 1))
155 radius_ *= rMultiplier;
156 radius_ =
max(radius_, (1.001*0.5*
mag(chord)));
167 const point newCentre =
180 calcFromCentre(p1, p3, newCentre,
false);
185 cs_ = coordSystem::cylindrical(centre, arcAxis, r1);
192 Foam::blockEdges::arcEdge::arcEdge
204 calcFromCentre(firstPoint(), lastPoint(), origin);
208 Foam::blockEdges::arcEdge::arcEdge
220 calcFromMidPoint(firstPoint(), lastPoint(),
midPoint);
224 Foam::blockEdges::arcEdge::arcEdge
236 Foam::blockEdges::arcEdge::arcEdge
248 Foam::blockEdges::arcEdge::arcEdge
271 scalar rMultiplier = 1;
276 rMultiplier = tok.number();
311 if (lambda < -SMALL || lambda > 1 + SMALL)
314 <<
"Limit parameter to [0-1] range: " <<
lambda <<
nl;
322 else if (
lambda >= 1 - SMALL)
327 return cs_.globalPosition(
vector(radius_, (
lambda*angle_), 0));
333 return (radius_*angle_);
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...
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.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Unit conversion functions.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
constexpr char nl
The newline '\n' character (0x0a)
A token holds an item read from Istream.
dimensionedScalar sqrt(const dimensionedScalar &ds)
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
virtual void write(Ostream &os) const =0
Write information.
Macros for easy insertion into run-time selection tables.
A blockEdge defined as an arc of a circle.
virtual const point & origin() const
Return origin.
void putBack(const token &tok)
Put back a token. Only a single put back is permitted.
point position(const scalar lambda) const
The point corresponding to the curve parameter [0-1].
Mid-point interpolation (weighting factors = 0.5) scheme class.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
constexpr scalar twoPi(2 *M_PI)
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
Container for searchableSurfaces. The collection is specified as a dictionary. For example...
constexpr scalar pi(M_PI)
virtual const coordinateRotation & rotation() const
The rotation specification.
const point & lastPoint() const
The location of the last point.
errorManip< error > abort(error &err)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
defineTypeNameAndDebug(arcEdge, 0)
addToRunTimeSelectionTable(blockEdge, arcEdge, Istream)
bool equal(const T &s1, const T &s2)
Compare two values for equality.
Define a curved edge that is parameterized for 0<lambda<1 between the start/end points.
vector point
Point is a vector.
const label start_
Index of the first point.
messageStream Info
Information stream (stdout output on master, null elsewhere)
const point & firstPoint() const
The location of the first point.
const label end_
Index of the last point.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
scalar length() const noexcept
The length of the curve.
#define InfoInFunction
Report an information message using Foam::Info.