37 const scalar
a = this->
a();
38 const scalar
b = this->
b();
39 const scalar
c = this->
c();
40 const scalar
d = this->
d();
44 <<
"Coefficients of the characteristic cubic polynomial:" <<
nl 60 const scalar
p = -(fma(-
a,
c, w) + fma(
b,
b/3.0, -w));
61 const scalar q =
b*
b*
b*scalar(2)/27 -
b*
c*
a/3 +
d*
a*
a;
62 const scalar numDiscr =
p*
p*
p/27 + q*q/4;
63 const scalar discr = (
mag(numDiscr) > VSMALL) ? numDiscr : 0;
66 const bool threeReal = discr < 0;
67 const bool oneRealTwoComplex = discr > 0;
69 (
mag(
p) >
sqrt(SMALL)) && !(threeReal || oneRealTwoComplex);
74 <<
"Numerical discriminant:" <<
tab << numDiscr <<
nl 75 <<
"Adjusted discriminant:" <<
tab << discr <<
nl 76 <<
"Number and types of the roots:" <<
nl 77 <<
"threeReal = " << threeReal <<
nl 78 <<
"oneRealTwoComplex = " << oneRealTwoComplex <<
nl 79 <<
"twoReal = " << twoReal <<
nl 80 <<
"oneReal = " << !(threeReal || oneRealTwoComplex || twoReal) <<
nl 84 static const scalar sqrt3 =
sqrt(3.0);
90 const scalar wCbRe = -q/2;
91 const scalar wCbIm =
sqrt(-discr);
92 const scalar wAbs =
cbrt(
hypot(wCbRe, wCbIm));
93 const scalar wArg =
atan2(wCbIm, wCbRe)/3;
94 const scalar wRe = wAbs*
cos(wArg);
95 const scalar wIm = wAbs*
sin(wArg);
99 x = -wRe -
mag(wIm)*sqrt3 -
b/3;
106 else if (oneRealTwoComplex)
108 const scalar wCb = -q/2 -
sign(q)*
sqrt(discr);
109 const scalar w =
cbrt(wCb);
110 const scalar t = w -
p/(3*w);
118 const scalar xRe = -t/2 -
b/3;
119 const scalar xIm = sqrt3/2*(w +
p/3/w);
155 <<
"#######" <<
endl;
dimensionedScalar sign(const dimensionedScalar &ds)
Roots< 2 > roots() const
Return the roots of the quadratic equation with no particular order.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
scalar c() const noexcept
constexpr char tab
The tab '\t' character(0x09)
Templated storage for the roots of polynomial equations, plus flags to indicate the nature of the roo...
dimensionedScalar cos(const dimensionedScalar &ds)
scalar b() const noexcept
dimensionedScalar cbrt(const dimensionedScalar &ds)
scalar a() const noexcept
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
Container to encapsulate various operations for quadratic equation of the forms with real coefficient...
Container to encapsulate various operations for linear equation of the forms with real coefficients: ...
Roots< 3 > roots() const
Return the roots of the cubic equation with no particular order.
void type(const direction i, const roots::type t)
Set the type of the i-th root.
Roots< 1 > roots() const
Return the real root of the linear equation.
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensionedScalar hypot(const dimensionedScalar &x, const dimensionedScalar &y)
scalar d() const noexcept