mirror of
https://github.com/UpsilonNumworks/Upsilon.git
synced 2026-03-22 07:10:40 +01:00
503 lines
30 KiB
C++
503 lines
30 KiB
C++
#include <poincare/trigonometry.h>
|
|
#include <poincare/constant.h>
|
|
#include <poincare/symbol.h>
|
|
#include <poincare/preferences.h>
|
|
#include <poincare/undefined.h>
|
|
#include <poincare/rational.h>
|
|
#include <poincare/multiplication.h>
|
|
#include <poincare/sign_function.h>
|
|
#include <poincare/subtraction.h>
|
|
#include <poincare/derivative.h>
|
|
#include <poincare/decimal.h>
|
|
#include <poincare/float.h>
|
|
#include <poincare/power.h>
|
|
#include <poincare/addition.h>
|
|
#include <ion.h>
|
|
#include <assert.h>
|
|
#include <cmath>
|
|
#include <float.h>
|
|
|
|
namespace Poincare {
|
|
|
|
float Trigonometry::characteristicXRange(const Expression & e, Context & context, Preferences::AngleUnit angleUnit) {
|
|
assert(e.numberOfChildren() == 1);
|
|
const char x[] = {Symbol::SpecialSymbols::UnknownX, 0};
|
|
int d = e.childAtIndex(0).polynomialDegree(context, x);
|
|
if (d < 0 || d > 1) {
|
|
// child(0) is not linear so we cannot easily find an interesting range
|
|
return e.childAtIndex(0).characteristicXRange(context, angleUnit);
|
|
}
|
|
// The expression e is x-independent
|
|
if (d == 0) {
|
|
return 0.0f;
|
|
}
|
|
// e has the form cos/sin/tan(ax+b) so it is periodic of period 2*Pi/a
|
|
assert(d == 1);
|
|
/* To compute a, the slope of the expression child(0), we compute the
|
|
* derivative of child(0) for any x value. */
|
|
Poincare::Derivative derivative = Poincare::Derivative::Builder(e.childAtIndex(0).clone(), Symbol(x, 1), Float<float>(1.0f));
|
|
float a = derivative.approximateToScalar<float>(context, angleUnit);
|
|
float pi = angleUnit == Preferences::AngleUnit::Radian ? M_PI : 180.0f;
|
|
return 2.0f*pi/std::fabs(a);
|
|
}
|
|
|
|
Expression Trigonometry::shallowReduceDirectFunction(Expression & e, Context& context, Preferences::AngleUnit angleUnit, ExpressionNode::ReductionTarget target) {
|
|
assert(e.type() == ExpressionNode::Type::Sine
|
|
|| e.type() == ExpressionNode::Type::Cosine
|
|
|| e.type() == ExpressionNode::Type::Tangent);
|
|
|
|
// Step 1. Try finding an easy standard calculation reduction
|
|
Expression lookup = Trigonometry::table(e.childAtIndex(0), e.type(), context, angleUnit, target);
|
|
if (!lookup.isUninitialized()) {
|
|
e.replaceWithInPlace(lookup);
|
|
return lookup;
|
|
}
|
|
|
|
// Step 2. Look for an expression of type "cos(arccos(x))", return x
|
|
ExpressionNode::Type correspondingType = e.type() == ExpressionNode::Type::Cosine ? ExpressionNode::Type::ArcCosine :
|
|
(e.type() == ExpressionNode::Type::Sine ? ExpressionNode::Type::ArcSine : ExpressionNode::Type::ArcTangent);
|
|
if (e.childAtIndex(0).type() == correspondingType) {
|
|
Expression result = e.childAtIndex(0).childAtIndex(0);
|
|
e.replaceWithInPlace(result);
|
|
return result;
|
|
}
|
|
|
|
// Step 3. Look for an expression of type "cos(arcsin(x))" or "sin(arccos(x)), return sqrt(1-x^2)
|
|
if ((e.type() == ExpressionNode::Type::Cosine && e.childAtIndex(0).type() == ExpressionNode::Type::ArcSine) || (e.type() == ExpressionNode::Type::Sine && e.childAtIndex(0).type() == ExpressionNode::Type::ArcCosine)) {
|
|
Expression sqrt =
|
|
Power(
|
|
Addition(
|
|
Rational(1),
|
|
Multiplication(
|
|
Rational(-1),
|
|
Power(e.childAtIndex(0).childAtIndex(0), Rational(2)) )
|
|
),
|
|
Rational(1,2)
|
|
);
|
|
// reduce x^2
|
|
sqrt.childAtIndex(0).childAtIndex(1).childAtIndex(1).shallowReduce(context, angleUnit, target);
|
|
// reduce -1*x^2
|
|
sqrt.childAtIndex(0).childAtIndex(1).shallowReduce(context, angleUnit, target);
|
|
// reduce 1-1*x^2
|
|
sqrt.childAtIndex(0).shallowReduce(context, angleUnit, target);
|
|
e.replaceWithInPlace(sqrt);
|
|
// reduce sqrt(1+(-1)*x^2)
|
|
return sqrt.shallowReduce(context, angleUnit, target);
|
|
}
|
|
|
|
// Step 4. Look for an expression of type "cos(arctan(x))" or "sin(arctan(x))"
|
|
// cos(arctan(x)) --> 1/sqrt(1+x^2)
|
|
// sin(arctan(x)) --> x/sqrt(1+x^2)
|
|
if ((e.type() == ExpressionNode::Type::Cosine || e.type() == ExpressionNode::Type::Sine) && e.childAtIndex(0).type() == ExpressionNode::Type::ArcTangent) {
|
|
Expression x = e.childAtIndex(0).childAtIndex(0);
|
|
// Build 1/sqrt(1+x^2)
|
|
Expression res =
|
|
Power(
|
|
Addition(
|
|
Rational(1),
|
|
Power(
|
|
e.type() == ExpressionNode::Type::Cosine ? x : x.clone(),
|
|
Rational(2))
|
|
),
|
|
Rational(-1,2)
|
|
);
|
|
|
|
// reduce x^2
|
|
res.childAtIndex(0).childAtIndex(1).shallowReduce(context, angleUnit, target);
|
|
// reduce 1+*x^2
|
|
res.childAtIndex(0).shallowReduce(context, angleUnit, target);
|
|
if (e.type() == ExpressionNode::Type::Sine) {
|
|
res = Multiplication(x, res);
|
|
// reduce (1+x^2)^(-1/2)
|
|
res.childAtIndex(0).shallowReduce(context, angleUnit, target);
|
|
}
|
|
e.replaceWithInPlace(res);
|
|
// reduce (1+x^2)^(-1/2) or x*(1+x^2)^(-1/2)
|
|
return res.shallowReduce(context, angleUnit, target);
|
|
}
|
|
|
|
// Step 3. Look for an expression of type "cos(-a)", return "+/-cos(a)"
|
|
Expression positiveArg = e.childAtIndex(0).makePositiveAnyNegativeNumeralFactor(context, angleUnit, target);
|
|
if (!positiveArg.isUninitialized()) {
|
|
// The argument was of form cos(-a)
|
|
if (e.type() == ExpressionNode::Type::Cosine) {
|
|
// cos(-a) = cos(a)
|
|
return e.shallowReduce(context, angleUnit, target);
|
|
} else {
|
|
// sin(-a) = -sin(a) or tan(-a) = -tan(a)
|
|
Multiplication m(Rational(-1));
|
|
e.replaceWithInPlace(m);
|
|
m.addChildAtIndexInPlace(e, 1, 1);
|
|
e.shallowReduce(context, angleUnit, target);
|
|
return m.shallowReduce(context, angleUnit, target);
|
|
}
|
|
}
|
|
|
|
/* Step 4. Look for an expression of type "cos(p/q * Pi)" in radians or
|
|
* "cos(p/q)" in degrees, put the argument in [0, Pi/2[ or [0, 90[ and
|
|
* multiply the cos/sin/tan by -1 if needed.
|
|
* We know thanks to Step 3 that p/q > 0. */
|
|
if ((angleUnit == Preferences::AngleUnit::Radian
|
|
&& e.childAtIndex(0).type() == ExpressionNode::Type::Multiplication
|
|
&& e.childAtIndex(0).numberOfChildren() == 2
|
|
&& e.childAtIndex(0).childAtIndex(1).type() == ExpressionNode::Type::Constant
|
|
&& e.childAtIndex(0).childAtIndex(1).convert<Constant>().isPi()
|
|
&& e.childAtIndex(0).childAtIndex(0).type() == ExpressionNode::Type::Rational)
|
|
|| (angleUnit == Preferences::AngleUnit::Degree
|
|
&& e.childAtIndex(0).type() == ExpressionNode::Type::Rational))
|
|
{
|
|
Rational r = angleUnit == Preferences::AngleUnit::Radian ? e.childAtIndex(0).childAtIndex(0).convert<Rational>() : e.childAtIndex(0).convert<Rational>();
|
|
/* Step 4.1. In radians:
|
|
* We first check if p/q * Pi is already in the right quadrant:
|
|
* p/q * Pi < Pi/2 => p/q < 2 => 2p < q */
|
|
Integer dividand = angleUnit == Preferences::AngleUnit::Radian ? Integer::Addition(r.unsignedIntegerNumerator(), r.unsignedIntegerNumerator()) : r.unsignedIntegerNumerator();
|
|
Integer divisor = angleUnit == Preferences::AngleUnit::Radian ? r.integerDenominator() : Integer::Multiplication(r.integerDenominator(), Integer(90));
|
|
if (divisor.isLowerThan(dividand)) {
|
|
/* Step 4.2. p/q * Pi is not in the wanted trigonometrical quadrant.
|
|
* We could subtract n*Pi to p/q with n an integer.
|
|
* Given p/q = (q'*q+r')/q, we have
|
|
* (p/q * Pi - q'*Pi) < Pi/2 => r'/q < 1/2 => 2*r'<q
|
|
* (q' is the theoretical n).*/
|
|
int unaryCoefficient = 1; // store 1 or -1 for the final result.
|
|
Integer piDivisor = angleUnit == Preferences::AngleUnit::Radian ? r.integerDenominator() : Integer::Multiplication(r.integerDenominator(), Integer(180));
|
|
IntegerDivision div = Integer::Division(r.unsignedIntegerNumerator(), piDivisor);
|
|
dividand = angleUnit == Preferences::AngleUnit::Radian ? Integer::Addition(div.remainder, div.remainder) : div.remainder;
|
|
if (divisor.isLowerThan(dividand)) {
|
|
/* Step 4.3. r'/q * Pi is not in the wanted trigonometrical quadrant,
|
|
* and because r'<q (as r' is the remainder of an euclidian division
|
|
* by q), we know that r'/q*Pi is in [Pi/2; Pi[.
|
|
* So we can take the new angle Pi - r'/q*Pi, which changes cosinus or
|
|
* tangent, but not sinus. The new rational is 1-r'/q = (q-r')/q. */
|
|
div.remainder = Integer::Subtraction(piDivisor, div.remainder);
|
|
if (e.type() == ExpressionNode::Type::Cosine || e.type() == ExpressionNode::Type::Tangent) {
|
|
unaryCoefficient *= -1;
|
|
}
|
|
}
|
|
if (div.remainder.isInfinity()) {
|
|
return e;
|
|
}
|
|
// Step 4.5. Build the new result.
|
|
Integer rDenominator = r.integerDenominator();
|
|
Expression newR = Rational(div.remainder, rDenominator);
|
|
Expression rationalParent = angleUnit == Preferences::AngleUnit::Radian ? e.childAtIndex(0) : e;
|
|
rationalParent.replaceChildAtIndexInPlace(0, newR);
|
|
newR.shallowReduce(context, angleUnit, target);
|
|
if (angleUnit == Preferences::AngleUnit::Radian) {
|
|
e.childAtIndex(0).shallowReduce(context, angleUnit, target);
|
|
}
|
|
if (Integer::Division(div.quotient, Integer(2)).remainder.isOne() && e.type() != ExpressionNode::Type::Tangent) {
|
|
/* Step 4.6. If we subtracted an odd number of Pi in 4.2, we need to
|
|
* multiply the result by -1 (because cos((2k+1)Pi + x) = -cos(x) */
|
|
unaryCoefficient *= -1;
|
|
}
|
|
Expression simplifiedCosine = e.shallowReduce(context, angleUnit, target); // recursive
|
|
Multiplication m = Multiplication(Rational(unaryCoefficient));
|
|
simplifiedCosine.replaceWithInPlace(m);
|
|
m.addChildAtIndexInPlace(simplifiedCosine, 1, 1);
|
|
return m.shallowReduce(context, angleUnit, target);
|
|
}
|
|
assert(r.sign() == ExpressionNode::Sign::Positive);
|
|
}
|
|
return e;
|
|
}
|
|
|
|
bool Trigonometry::ExpressionIsEquivalentToTangent(const Expression & e) {
|
|
// We look for (cos^-1 * sin)
|
|
assert(ExpressionNode::Type::Power < ExpressionNode::Type::Sine);
|
|
if (e.type() == ExpressionNode::Type::Multiplication
|
|
&& e.childAtIndex(1).type() == ExpressionNode::Type::Sine
|
|
&& e.childAtIndex(0).type() == ExpressionNode::Type::Power
|
|
&& e.childAtIndex(0).childAtIndex(0).type() == ExpressionNode::Type::Cosine
|
|
&& e.childAtIndex(0).childAtIndex(1).type() == ExpressionNode::Type::Rational
|
|
&& e.childAtIndex(0).childAtIndex(1).convert<Rational>().isMinusOne()) {
|
|
return true;
|
|
}
|
|
return false;
|
|
}
|
|
|
|
bool Trigonometry::parentIsDirectTrigonometry(const Expression & e) {
|
|
Expression parent = e.parent();
|
|
if (parent.isUninitialized()) {
|
|
return false;
|
|
}
|
|
if (parent.type() == ExpressionNode::Type::Cosine || parent.type() == ExpressionNode::Type::Sine || parent.type() == ExpressionNode::Type::Tangent) {
|
|
return true;
|
|
}
|
|
return false;
|
|
}
|
|
|
|
Expression Trigonometry::shallowReduceInverseFunction(Expression & e, Context& context, Preferences::AngleUnit angleUnit, ExpressionNode::ReductionTarget target) {
|
|
assert(e.type() == ExpressionNode::Type::ArcCosine || e.type() == ExpressionNode::Type::ArcSine || e.type() == ExpressionNode::Type::ArcTangent);
|
|
ExpressionNode::Type correspondingType = e.type() == ExpressionNode::Type::ArcCosine ? ExpressionNode::Type::Cosine : (e.type() == ExpressionNode::Type::ArcSine ? ExpressionNode::Type::Sine : ExpressionNode::Type::Tangent);
|
|
float pi = angleUnit == Preferences::AngleUnit::Radian ? M_PI : 180;
|
|
|
|
// Step 1. Look for an expression of type "arccos(cos(x))", return x
|
|
if (e.childAtIndex(0).type() == correspondingType) {
|
|
float trigoOp = e.childAtIndex(0).childAtIndex(0).approximateToScalar<float>(context, angleUnit);
|
|
if ((e.type() == ExpressionNode::Type::ArcCosine && trigoOp >= 0.0f && trigoOp <= pi) ||
|
|
(e.type() == ExpressionNode::Type::ArcSine && trigoOp >= -pi/2.0f && trigoOp <= pi/2.0f) ||
|
|
(e.type() == ExpressionNode::Type::ArcTangent && trigoOp >= -pi/2.0f && trigoOp <= pi/2.0f)) {
|
|
Expression result = e.childAtIndex(0).childAtIndex(0);
|
|
e.replaceWithInPlace(result);
|
|
return result;
|
|
}
|
|
}
|
|
|
|
// Step 2. Special case for arctan(sin(x)/cos(x))
|
|
if (e.type() == ExpressionNode::Type::ArcTangent && ExpressionIsEquivalentToTangent(e.childAtIndex(0))) {
|
|
float trigoOp = e.childAtIndex(0).childAtIndex(1).childAtIndex(0).approximateToScalar<float>(context, angleUnit);
|
|
if (trigoOp >= -pi/2.0f && trigoOp <= pi/2.0f) {
|
|
Expression result = e.childAtIndex(0).childAtIndex(1).childAtIndex(0);
|
|
e.replaceWithInPlace(result);
|
|
return result;
|
|
}
|
|
}
|
|
|
|
// Step 3. Look for an expression of type "arctan(1/X), return sign(x)*Pi/2-arctan(x)
|
|
if (e.type() == ExpressionNode::Type::ArcTangent && e.childAtIndex(0).type() == ExpressionNode::Type::Power && e.childAtIndex(0).childAtIndex(1).type() == ExpressionNode::Type::Rational && e.childAtIndex(0).childAtIndex(1).convert<Rational>().isMinusOne()) {
|
|
Expression x = e.childAtIndex(0).childAtIndex(0);
|
|
Expression sign = SignFunction::Builder(x.clone());
|
|
Multiplication m0(Rational(1,2), sign, Constant(Ion::Charset::SmallPi));
|
|
sign.shallowReduce(context, angleUnit, target);
|
|
e.replaceChildAtIndexInPlace(0, x);
|
|
Addition a(m0);
|
|
e.replaceWithInPlace(a);
|
|
Multiplication m1(Rational(-1), e);
|
|
e.shallowReduce(context, angleUnit, target);
|
|
a.addChildAtIndexInPlace(m1, 1, 1);
|
|
return a.shallowReduce(context, angleUnit, target);
|
|
}
|
|
|
|
// Step 4. Try finding an easy standard calculation reduction
|
|
Expression lookup = Trigonometry::table(e.childAtIndex(0), e.type(), context, angleUnit, target);
|
|
if (!lookup.isUninitialized()) {
|
|
e.replaceWithInPlace(lookup);
|
|
return lookup;
|
|
}
|
|
|
|
/* We do not apply some rules if:
|
|
* - the parent node is a cosine, a sine or a tangent. In this case there is a simplication of
|
|
* form f(g(x)) with f cos, sin or tan and g acos, asin or atan.
|
|
* - the reduction is being BottomUp. In this case, we do not yet have any
|
|
* information on the parent which could later be a cosine, a sine or a tangent.
|
|
*/
|
|
bool letArcFunctionAtRoot = target == ExpressionNode::ReductionTarget::BottomUpComputation || parentIsDirectTrigonometry(e);
|
|
/* Step 5. Handle opposite argument: arccos(-x) = Pi-arcos(x),
|
|
* arcsin(-x) = -arcsin(x), arctan(-x)= -arctan(x) *
|
|
*/
|
|
if (!letArcFunctionAtRoot) {
|
|
Expression positiveArg = e.childAtIndex(0).makePositiveAnyNegativeNumeralFactor(context, angleUnit, target);
|
|
if (!positiveArg.isUninitialized()) {
|
|
// The argument was made positive
|
|
// acos(-x) = pi-acos(x)
|
|
if (e.type() == ExpressionNode::Type::ArcCosine) {
|
|
Expression pi = angleUnit == Preferences::AngleUnit::Radian ? static_cast<Expression>(Constant(Ion::Charset::SmallPi)) : static_cast<Expression>(Rational(180));
|
|
Subtraction s;
|
|
e.replaceWithInPlace(s);
|
|
s.replaceChildAtIndexInPlace(0, pi);
|
|
s.replaceChildAtIndexInPlace(1, e);
|
|
e.shallowReduce(context, angleUnit, target);
|
|
return s.shallowReduce(context, angleUnit, target);
|
|
} else {
|
|
// asin(-x) = -asin(x) or atan(-x) = -atan(x)
|
|
Multiplication m(Rational(-1));
|
|
e.replaceWithInPlace(m);
|
|
m.addChildAtIndexInPlace(e, 1, 1);
|
|
e.shallowReduce(context, angleUnit, target);
|
|
return m.shallowReduce(context, angleUnit, target);
|
|
}
|
|
}
|
|
}
|
|
|
|
return e;
|
|
}
|
|
|
|
/* We use the cheat table to look for known simplifications (e.g.
|
|
* cos(0)=1, cos(Pi/2)=1...). For each entry of the table, we store its
|
|
* expression and its float approximation in order to quickly scan the table
|
|
* looking for our input approximation. If one entry matches, we then check that
|
|
* the actual expression of our input is equivalent to the table expression.
|
|
*/
|
|
|
|
static_assert('\x8A' == Ion::Charset::SmallPi, "Unicode error");
|
|
|
|
struct Pair {
|
|
constexpr Pair(const char * e, float v = NAN) : expression(e), value(v) {}
|
|
const char * expression;
|
|
float value;
|
|
};
|
|
constexpr Pair cheatTable[Trigonometry::k_numberOfEntries][5] =
|
|
// Angle in Radian | Angle in Degree | Cosine | Sine | Tangent
|
|
{{Pair("-90", -90.0f), Pair("\x8A*(-2)^(-1)", -1.5707963267948966f), "",
|
|
Pair("-1",-1.0f), "undef"},
|
|
{Pair("-75",-75.0), Pair("\x8A*(-5)*12^(-1)",-1.3089969389957472f), "",
|
|
Pair("(-1)*6^(1/2)*4^(-1)-2^(1/2)*4^(-1)",-0.9659258262890683f), Pair("-(3^(1/2)+2)",-3.7320508075688776f)},
|
|
{Pair("-72",-72.0), Pair("\x8A*2*(-5)^(-1)",-1.2566370614359172f), "",
|
|
Pair("-(5/8+5^(1/2)/8)^(1/2)",-0.9510565162951535f), Pair("-(5+2*5^(1/2))^(1/2)",-3.077683537175253f)},
|
|
{Pair("-135/2",67.5f), Pair("\x8A*(-3)*8^(-1)",-1.1780972450961724f), "",
|
|
Pair("-(2+2^(1/2))^(1/2)*2^(-1)",-0.9238795325112867f), Pair("-1-2^(1/2)",-2.4142135623730945f)},
|
|
{Pair("-60",-60.0f), Pair("\x8A*(-3)^(-1)",-1.0471975511965976f), "",
|
|
Pair("-3^(1/2)*2^(-1)",-0.8660254037844386f), Pair("-3^(1/2)",-1.7320508075688767f)},
|
|
{Pair("-54",-54.0f), Pair("\x8A*(-3)*10^(-1)",-0.9424777960769379), "",
|
|
Pair("4^(-1)*(-1-5^(1/2))",-0.8090169943749473f), Pair("-(1+2*5^(-1/2))^(1/2)",-1.3763819204711731f)},
|
|
{Pair("-45",-45.0f), Pair("\x8A*(-4)^(-1)",-0.7853981633974483f), "",
|
|
Pair("(-1)*(2^(-1/2))",-0.7071067811865475f), Pair("-1",-1.0f)},
|
|
{Pair("-36",-36.0f), Pair("\x8A*(-5)^(-1)",-0.6283185307179586f), "",
|
|
Pair("-(5/8-5^(1/2)/8)^(1/2)",-0.5877852522924731f), Pair("-(5-2*5^(1/2))^(1/2)",-0.7265425280053609f)},
|
|
{Pair("-30",-30.0f), Pair("\x8A*(-6)^(-1)",-0.5235987755982988f), "",
|
|
Pair("-0.5",-0.5f), Pair("-3^(-1/2)",-0.5773502691896256f)},
|
|
{Pair("-45/2",-22.5f), Pair("\x8A*(-8)^(-1)",-0.39269908169872414f), "",
|
|
Pair("(2-2^(1/2))^(1/2)*(-2)^(-1)",-0.3826834323650898f), Pair("1-2^(1/2)",-0.4142135623730951f)},
|
|
{Pair("-18",-18.0f), Pair("\x8A*(-10)^(-1)",-0.3141592653589793f), "",
|
|
Pair("4^(-1)*(1-5^(1/2))",-0.3090169943749474f), Pair("-(1-2*5^(-1/2))^(1/2)",-0.3249196962329063f)},
|
|
{Pair("-15",-15.0f), Pair("\x8A*(-12)^(-1)",-0.2617993877991494f), "",
|
|
Pair("-6^(1/2)*4^(-1)+2^(1/2)*4^(-1)",-0.25881904510252074f), Pair("3^(1/2)-2",-0.2679491924311227f)},
|
|
{Pair("0",0.0f), Pair("0",0.0f), Pair("1",1.0f),
|
|
Pair("0",0.0f), Pair("0",0.0f)},
|
|
{Pair("15",15.0f), Pair("\x8A*12^(-1)",0.2617993877991494f), Pair("6^(1/2)*4^(-1)+2^(1/2)*4^(-1)",0.9659258262890683f),
|
|
Pair("6^(1/2)*4^(-1)+2^(1/2)*(-4)^(-1)",0.25881904510252074f), Pair("-(3^(1/2)-2)",0.2679491924311227f)},
|
|
{Pair("18",18.0f), Pair("\x8A*10^(-1)",0.3141592653589793f), Pair("(5/8+5^(1/2)/8)^(1/2)",0.9510565162951535f),
|
|
Pair("4^(-1)*(5^(1/2)-1)",0.3090169943749474f), Pair("(1-2*5^(-1/2))^(1/2)",0.3249196962329063f)},
|
|
{Pair("45/2",22.5f), Pair("\x8A*8^(-1)",0.39269908169872414f), Pair("(2+2^(1/2))^(1/2)*2^(-1)",0.9238795325112867f),
|
|
Pair("(2-2^(1/2))^(1/2)*2^(-1)",0.3826834323650898f), Pair("2^(1/2)-1",0.4142135623730951f)},
|
|
{Pair("30",30.0f), Pair("\x8A*6^(-1)",0.5235987755982988f), Pair("3^(1/2)*2^(-1)",0.8660254037844387f),
|
|
Pair("0.5",0.5f), Pair("3^(-1/2)",0.5773502691896256f)},
|
|
{Pair("36",36.0f), Pair("\x8A*5^(-1)",0.6283185307179586f), Pair("(5^(1/2)+1)*4^(-1)",0.8090169943749475f),
|
|
Pair("(5/8-5^(1/2)/8)^(1/2)",0.5877852522924731f), Pair("(5-2*5^(1/2))^(1/2)",0.7265425280053609f)},
|
|
{Pair("45",45.0f), Pair("\x8A*4^(-1)",0.7853981633974483f), Pair("2^(-1/2)",0.7071067811865476f),
|
|
Pair("2^(-1/2)",0.7071067811865475f), Pair("1",1.0f)},
|
|
{Pair("54",54.0f), Pair("\x8A*3*10^(-1)",0.9424777960769379f), Pair("(5/8-5^(1/2)/8)^(1/2)",0.5877852522924732f),
|
|
Pair("4^(-1)*(5^(1/2)+1)",0.8090169943749473f), Pair("(1+2*5^(-1/2))^(1/2)",1.3763819204711731f)},
|
|
{Pair("60",60.0f), Pair("\x8A*3^(-1)",1.0471975511965976f), Pair("0.5",0.5f),
|
|
Pair("3^(1/2)*2^(-1)",0.8660254037844386f), Pair("3^(1/2)",1.7320508075688767f)},
|
|
{Pair("135/2",67.5f), Pair("\x8A*3*8^(-1)",1.1780972450961724f), Pair("(2-2^(1/2))^(1/2)*2^(-1)",0.38268343236508984f),
|
|
Pair("(2+2^(1/2))^(1/2)*2^(-1)",0.9238795325112867f), Pair("1+2^(1/2)",2.4142135623730945f)},
|
|
{Pair("72",72.0f), Pair("\x8A*2*5^(-1)",1.2566370614359172f), Pair("(5^(1/2)-1)*4^(-1)",0.30901699437494745f),
|
|
Pair("(5/8+5^(1/2)/8)^(1/2)",0.9510565162951535f), Pair("(5+2*5^(1/2))^(1/2)",3.077683537175253f)},
|
|
{Pair("75",75.0f), Pair("\x8A*5*12^(-1)",1.3089969389957472f), Pair("6^(1/2)*4^(-1)+2^(1/2)*(-4)^(-1)",0.25881904510252074f),
|
|
Pair("6^(1/2)*4^(-1)+2^(1/2)*4^(-1)",0.9659258262890683f), Pair("3^(1/2)+2",3.7320508075688776f)},
|
|
{Pair("90",90.0f), Pair("\x8A*2^(-1)",1.5707963267948966f), Pair("0",0.0f),
|
|
Pair("1",1.0f), "undef"},
|
|
{Pair("105",105.0f), Pair("\x8A*7*12^(-1)",1.832595714594046f), Pair("-6^(1/2)*4^(-1)+2^(1/2)*4^(-1)",-0.25881904510252063f),
|
|
"", ""},
|
|
{Pair("108",108.0f), Pair("\x8A*3*5^(-1)",1.8849555921538759f), Pair("(1-5^(1/2))*4^(-1)",-0.30901699437494734f),
|
|
"", ""},
|
|
{Pair("225/2",112.5f), Pair("\x8A*5*8^(-1)",1.9634954084936207f), Pair("(2-2^(1/2))^(1/2)*(-2)^(-1)",-0.3826834323650897f),
|
|
"", ""},
|
|
{Pair("120",120.0f), Pair("\x8A*2*3^(-1)",2.0943951023931953f), Pair("-0.5",-0.5f),
|
|
"", ""},
|
|
{Pair("126",126.0f), Pair("\x8A*7*10^(-1)",2.199114857512855f), Pair("-(5*8^(-1)-5^(1/2)*8^(-1))^(1/2)",-0.587785252292473f),
|
|
"", ""},
|
|
{Pair("135",135.0f), Pair("\x8A*3*4^(-1)",2.356194490192345f), Pair("(-1)*(2^(-1/2))",-0.7071067811865475f),
|
|
"", ""},
|
|
{Pair("144",144.0f), Pair("\x8A*4*5^(-1)",2.5132741228718345f), Pair("(-5^(1/2)-1)*4^(-1)",-0.8090169943749473f),
|
|
"", ""},
|
|
{Pair("150",150.0f), Pair("\x8A*5*6^(-1)",2.6179938779914944f), Pair("-3^(1/2)*2^(-1)",-0.8660254037844387f),
|
|
"", ""},
|
|
{Pair("315/2",157.5f), Pair("\x8A*7*8^(-1)",2.748893571891069f), Pair("-(2+2^(1/2))^(1/2)*2^(-1)",-0.9238795325112867f),
|
|
"", ""},
|
|
{Pair("162",162.0f), Pair("\x8A*9*10^(-1)",2.827433388230814f), Pair("-(5*8^(-1)+5^(1/2)*8^(-1))^(1/2)",-0.9510565162951535f),
|
|
"", ""},
|
|
{Pair("165",165.0f), Pair("\x8A*11*12^(-1)",2.8797932657906435f), Pair("(-1)*6^(1/2)*4^(-1)-2^(1/2)*4^(-1)",-0.9659258262890682f),
|
|
"", ""},
|
|
{Pair("180",180.0f), Pair("\x8A",3.141592653589793f), Pair("-1",-1.0f),
|
|
Pair("0",0.0f), Pair("0",0.0f)}};
|
|
|
|
Expression Trigonometry::table(const Expression e, ExpressionNode::Type type, Context & context, Preferences::AngleUnit angleUnit, ExpressionNode::ReductionTarget target) {
|
|
assert(type == ExpressionNode::Type::Sine
|
|
|| type == ExpressionNode::Type::Cosine
|
|
|| type == ExpressionNode::Type::Tangent
|
|
|| type == ExpressionNode::Type::ArcCosine
|
|
|| type == ExpressionNode::Type::ArcSine
|
|
|| type == ExpressionNode::Type::ArcTangent);
|
|
|
|
int angleUnitIndex = angleUnit == Preferences::AngleUnit::Radian ? 1 : 0;
|
|
int trigonometricFunctionIndex = type == ExpressionNode::Type::Cosine || type == ExpressionNode::Type::ArcCosine ? 2 : (type == ExpressionNode::Type::Sine || type == ExpressionNode::Type::ArcSine ? 3 : 4);
|
|
int inputIndex = type == ExpressionNode::Type::ArcCosine || type == ExpressionNode::Type::ArcSine || type == ExpressionNode::Type::ArcTangent ? trigonometricFunctionIndex : angleUnitIndex;
|
|
int outputIndex = type == ExpressionNode::Type::ArcCosine || type == ExpressionNode::Type::ArcSine || type == ExpressionNode::Type::ArcTangent ? angleUnitIndex : trigonometricFunctionIndex;
|
|
|
|
/* Avoid looping if we can exclude quickly that the e is in the table */
|
|
if (inputIndex == 0 && e.type() != ExpressionNode::Type::Rational) {
|
|
return Expression();
|
|
}
|
|
if (inputIndex == 1 && e.type() != ExpressionNode::Type::Rational && e.type() != ExpressionNode::Type::Multiplication && e.type() != ExpressionNode::Type::Constant) {
|
|
return Expression();
|
|
}
|
|
if (inputIndex >1 && e.type() != ExpressionNode::Type::Rational && e.type() != ExpressionNode::Type::Multiplication && e.type() != ExpressionNode::Type::Power && e.type() != ExpressionNode::Type::Addition) {
|
|
return Expression();
|
|
}
|
|
// We approximate the given expression to quickly compare it to the cheat table entries.
|
|
float eValue = e.approximateToScalar<float>(context, angleUnit);
|
|
if (std::isnan(eValue) || std::isinf(eValue)) {
|
|
return Expression();
|
|
}
|
|
for (int i = 0; i < k_numberOfEntries; i++) {
|
|
float inputValue = cheatTable[i][inputIndex].value;
|
|
if (std::isnan(inputValue) || std::fabs(inputValue-eValue) > Expression::epsilon<float>()) {
|
|
continue;
|
|
}
|
|
// Our given expression approximation matches a table entry, we check that both expressions are identical
|
|
Expression input = Expression::Parse(cheatTable[i][inputIndex].expression);
|
|
assert(!input.isUninitialized());
|
|
input = input.deepReduce(context, angleUnit, target);
|
|
bool rightInput = input.isIdenticalTo(e);
|
|
if (rightInput) {
|
|
Expression output = Expression::Parse(cheatTable[i][outputIndex].expression);
|
|
if (output.isUninitialized()) {
|
|
return Expression();
|
|
}
|
|
return output.deepReduce(context, angleUnit, target);
|
|
}
|
|
break;
|
|
}
|
|
return Expression();
|
|
}
|
|
|
|
template <typename T>
|
|
std::complex<T> Trigonometry::ConvertToRadian(const std::complex<T> c, Preferences::AngleUnit angleUnit) {
|
|
if (angleUnit == Preferences::AngleUnit::Degree) {
|
|
return c * std::complex<T>(M_PI/180.0);
|
|
}
|
|
return c;
|
|
}
|
|
|
|
template <typename T>
|
|
std::complex<T> Trigonometry::ConvertRadianToAngleUnit(const std::complex<T> c, Preferences::AngleUnit angleUnit) {
|
|
if (angleUnit == Preferences::AngleUnit::Degree) {
|
|
return c * std::complex<T>(180/M_PI);
|
|
}
|
|
return c;
|
|
}
|
|
|
|
template<typename T>
|
|
T Trigonometry::RoundToMeaningfulDigits(T result, T input) {
|
|
/* Cheat: openbsd trigonometric functions are numerical implementation and
|
|
* thus are approximative.
|
|
* The error epsilon is ~1E-7 on float and ~1E-15 on double. In order to
|
|
* avoid weird results as acos(1) = 6E-17 or cos(Pi/2) = 4E-17, we round
|
|
* the result to its 1E-6 or 1E-14 precision when its ratio with the
|
|
* argument (pi/2 in the exemple) is smaller than epsilon. This way, we
|
|
* have sin(pi) ~ 0 and sin(1E-15)=1E-15.
|
|
* We can't do that for all evaluation as the user can operate on values as
|
|
* small as 1E-308 (in double) and most results still be correct. */
|
|
if (input == 0.0 || std::fabs(result/input) <= Expression::epsilon<T>()) {
|
|
T precision = 10*Expression::epsilon<T>();
|
|
result = std::round(result/precision)*precision;
|
|
}
|
|
return result;
|
|
}
|
|
|
|
template <typename T>
|
|
std::complex<T> Trigonometry::RoundToMeaningfulDigits(const std::complex<T> result, const std::complex<T> input) {
|
|
return std::complex<T>(RoundToMeaningfulDigits(result.real(), std::abs(input)), RoundToMeaningfulDigits(result.imag(), std::abs(input)));
|
|
}
|
|
|
|
template std::complex<float> Trigonometry::ConvertToRadian<float>(std::complex<float>, Preferences::AngleUnit);
|
|
template std::complex<double> Trigonometry::ConvertToRadian<double>(std::complex<double>, Preferences::AngleUnit);
|
|
template std::complex<float> Trigonometry::ConvertRadianToAngleUnit<float>(std::complex<float>, Preferences::AngleUnit);
|
|
template std::complex<double> Trigonometry::ConvertRadianToAngleUnit<double>(std::complex<double>, Preferences::AngleUnit);
|
|
template std::complex<float> Trigonometry::RoundToMeaningfulDigits<float>(std::complex<float>, std::complex<float>);
|
|
template std::complex<double> Trigonometry::RoundToMeaningfulDigits<double>(std::complex<double>, std::complex<double>);
|
|
|
|
}
|