首页 > 代码库 > 一元二次方程求解C++实现
一元二次方程求解C++实现
-----------------------------------------------------------------------------------------------------------------------------
typedef double Number
int CesiumMath::sign(Number value){
if (value > 0) {
return 1;
}
if (value < 0) {
return -1;
}
return 0;
}
------------------------------------------------------------------------------------------------------------------------------
class QuadraticRealPolynomial{
public:
static Number computeDiscriminant(Number a, Number b, Number c);
static std::vector<Number>* computeRealRoots(Number a, Number b, Number c, std::vector<Number>& roots);
private:
static Number addWithCancellationCheck(Number left, Number right, Number tolerance);
public:
Number _root1;
Number _root2;
------------------------------------------------------------------------------------------------------------------------------
Number QuadraticRealPolynomial::computeDiscriminant(Number a, Number b, Number c)
{Number discriminant = b * b - 4.0 * a * c;
return discriminant;
}
std::vector<Number>* QuadraticRealPolynomial::computeRealRoots(Number a, Number b, Number c, std::vector<Number>& roots)
{
Number root1, root2;
roots.clear();
Number ratio;
if (a == 0.0)
{
if (b == 0.0)
{
// Constant function: c = 0.
return &roots;
}
// Linear function: b * x + c = 0.
root1 = -c/b;
roots.push_back(root1);
return &roots;
}
else if (b == 0.0)
{
if (c == 0.0)
{
// 2nd order monomial: a * x^2 = 0.
root1 = 0.0;
root2 = 0.0;
roots.push_back(root1);
roots.push_back(root2);
return &roots;
}
Number cMagnitude = std::abs(c);
Number aMagnitude = std::abs(a);
if ((cMagnitude < aMagnitude) && (cMagnitude / aMagnitude < CesiumMath::_EPSILON14))
{
// c ~= 0.0.
// 2nd order monomial: a * x^2 = 0.
root1 = 0.0;
root2 = 0.0;
roots.push_back(root1);
roots.push_back(root2);
return &roots;
} else if ((cMagnitude > aMagnitude) && (aMagnitude / cMagnitude < CesiumMath::_EPSILON14))
{
// a ~= 0.0.
// Constant function: c = 0.
return &roots;
}
// a * x^2 + c = 0
ratio = -c / a;
if (ratio < 0.0)
{
// Both roots are complex.
return &roots;
}
// Both roots are real.
Number root = std::sqrt(ratio);
root1 = -root;
root2 = root;
roots.push_back(root1);
roots.push_back(root2);
return &roots;
}
else if (c == 0.0)
{
// a * x^2 + b * x = 0
ratio = -b / a;
if (ratio < 0.0)
{
root1 = ratio;
root2 = 0.0;
roots.push_back(root1);
roots.push_back(root2);
return &roots;
}
root1 = 0.0;
root2 = ratio;
roots.push_back(root1);
roots.push_back(root2);
return &roots;
}
// a * x^2 + b * x + c = 0
Number b2 = b * b;
Number four_ac = 4.0 * a * c;
Number radicand = addWithCancellationCheck(b2, -four_ac, CesiumMath::_EPSILON14);
if (radicand < 0.0)
{
// Both roots are complex.
return &roots;
}
Number q = -0.5 * addWithCancellationCheck(b, CesiumMath::sign(b) * std::sqrt(radicand), CesiumMath::_EPSILON14);
if (b > 0.0)
{
root1 = q / a;
root2 = c / q;
roots.push_back(root1);
roots.push_back(root2);
return &roots;
}
root1 = c / q;
root2 = q / a;
roots.push_back(root1);
roots.push_back(root2);
return &roots;
}
Number QuadraticRealPolynomial::addWithCancellationCheck(Number left, Number right, Number tolerance)
{
Number difference = left + right;
if ((CesiumMath::sign(left) != CesiumMath::sign(right)) &&
std::abs(difference / std::max(std::abs(left), std::abs(right))) < tolerance)
{
return 0.0;
}
return difference;
}
一元二次方程求解C++实现