From 0316702adf3ac26d4dffc09fcb63e20b88d69d80 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9a=20Saviot?= Date: Fri, 23 Aug 2019 14:14:44 +0200 Subject: [PATCH] [poincare] NormCDF2 --- poincare/Makefile | 1 + poincare/include/poincare/expression.h | 6 ++ poincare/include/poincare/expression_node.h | 1 + poincare/include/poincare/norm_cdf2.h | 52 ++++++++++++++++ .../include/poincare/normal_distribution.h | 4 +- poincare/include/poincare_nodes.h | 1 + poincare/src/norm_cdf2.cpp | 62 +++++++++++++++++++ poincare/src/normal_distribution.cpp | 14 +++-- poincare/src/parsing/parser.h | 1 + poincare/src/tree_handle.cpp | 1 + poincare/test/approximation.cpp | 3 + 11 files changed, 139 insertions(+), 7 deletions(-) create mode 100644 poincare/include/poincare/norm_cdf2.h create mode 100644 poincare/src/norm_cdf2.cpp diff --git a/poincare/Makefile b/poincare/Makefile index 345d01750..ea7a78291 100644 --- a/poincare/Makefile +++ b/poincare/Makefile @@ -96,6 +96,7 @@ poincare_src += $(addprefix poincare/src/,\ n_ary_expression.cpp \ naperian_logarithm.cpp \ norm_cdf.cpp \ + norm_cdf2.cpp \ nth_root.cpp \ number.cpp \ opposite.cpp \ diff --git a/poincare/include/poincare/expression.h b/poincare/include/poincare/expression.h index 56192e156..5d4ef4662 100644 --- a/poincare/include/poincare/expression.h +++ b/poincare/include/poincare/expression.h @@ -66,6 +66,7 @@ class Expression : public TreeHandle { friend class MultiplicationNode; friend class NaperianLogarithm; friend class NormCDF; + friend class NormCDF2; friend class NthRoot; friend class Number; friend class Opposite; @@ -280,6 +281,11 @@ protected: assert(children.type() == ExpressionNode::Type::Matrix); return U::Builder(children.childAtIndex(0), children.childAtIndex(1), children.childAtIndex(2)); } + template + static Expression UntypedBuilderFourChildren(Expression children) { + assert(children.type() == ExpressionNode::Type::Matrix); + return U::Builder(children.childAtIndex(0), children.childAtIndex(1), children.childAtIndex(2), children.childAtIndex(3)); + } template T convert() const { /* This function allows to convert Expression to derived Expressions. diff --git a/poincare/include/poincare/expression_node.h b/poincare/include/poincare/expression_node.h index 15ce80ed6..c8b787151 100644 --- a/poincare/include/poincare/expression_node.h +++ b/poincare/include/poincare/expression_node.h @@ -73,6 +73,7 @@ public: MatrixTrace, NaperianLogarithm, NormCDF, + NormCDF2, NthRoot, Opposite, Parenthesis, diff --git a/poincare/include/poincare/norm_cdf2.h b/poincare/include/poincare/norm_cdf2.h new file mode 100644 index 000000000..a6e2d183d --- /dev/null +++ b/poincare/include/poincare/norm_cdf2.h @@ -0,0 +1,52 @@ +#ifndef POINCARE_NORMCDF2_H +#define POINCARE_NORMCDF2_H + +#include +#include + +namespace Poincare { + +class NormCDF2Node final : public ExpressionNode { +public: + + // TreeNode + size_t size() const override { return sizeof(NormCDF2Node); } + int numberOfChildren() const override; +#if POINCARE_TREE_LOG + virtual void logNodeName(std::ostream & stream) const override { + stream << "NormCDF2"; + } +#endif + + // Properties + Type type() const override { return Type::NormCDF2; } + Sign sign(Context * context) const override { return Sign::Positive; } + Expression setSign(Sign s, ReductionContext reductionContext) override; + +private: + // Layout + Layout createLayout(Preferences::PrintFloatMode floatDisplayMode, int numberOfSignificantDigits) const override; + int serialize(char * buffer, int bufferSize, Preferences::PrintFloatMode floatDisplayMode, int numberOfSignificantDigits) const override; + + // Simplication + Expression shallowReduce(ReductionContext reductionContext) override; + LayoutShape leftLayoutShape() const override { return LayoutShape::MoreLetters; }; + LayoutShape rightLayoutShape() const override { return LayoutShape::BoundaryPunctuation; } + + // Evaluation + Evaluation approximate(SinglePrecision p, Context * context, Preferences::ComplexFormat complexFormat, Preferences::AngleUnit angleUnit) const override { return templatedApproximate(context, complexFormat, angleUnit); } + Evaluation approximate(DoublePrecision p, Context * context, Preferences::ComplexFormat complexFormat, Preferences::AngleUnit angleUnit) const override { return templatedApproximate(context, complexFormat, angleUnit); } + template Evaluation templatedApproximate(Context * context, Preferences::ComplexFormat complexFormat, Preferences::AngleUnit angleUnit) const; +}; + +class NormCDF2 final : public Expression { +public: + NormCDF2(const NormCDF2Node * n) : Expression(n) {} + static NormCDF2 Builder(Expression child0, Expression child1, Expression child2, Expression child3) { return TreeHandle::FixedArityBuilder(ArrayBuilder(child0, child1, child2, child3).array(), 4); } + static constexpr Expression::FunctionHelper s_functionHelper = Expression::FunctionHelper("normcdf2", 4, &UntypedBuilderFourChildren); + Expression shallowReduce(ExpressionNode::ReductionContext reductionContext); +}; + +} + +#endif diff --git a/poincare/include/poincare/normal_distribution.h b/poincare/include/poincare/normal_distribution.h index bf4517be1..64f8d60c7 100644 --- a/poincare/include/poincare/normal_distribution.h +++ b/poincare/include/poincare/normal_distribution.h @@ -9,7 +9,7 @@ class NormalDistribution final { public: template static T EvaluateAtAbscissa(T x, T mu, T sigma); template static T CumulativeDistributiveFunctionAtAbscissa(T x, T mu, T sigma); - static double CumulativeDistributiveInverseForProbability(double probability, float mu, float sigma); + template static T CumulativeDistributiveInverseForProbability(T probability, T mu, T sigma); private: /* For the standard normal distribution, P(X < y) > 0.9999995 for y >= 4.892 so the * value displayed is 1. But this is dependent on the fact that we display @@ -17,7 +17,7 @@ private: static_assert(Preferences::LargeNumberOfSignificantDigits == 7, "k_boundStandardNormalDistribution is ill-defined compared to LargeNumberOfSignificantDigits"); constexpr static double k_boundStandardNormalDistribution = 4.892; template static T StandardNormalCumulativeDistributiveFunctionAtAbscissa(T abscissa); - static double StandardNormalCumulativeDistributiveInverseForProbability(double probability); + template static T StandardNormalCumulativeDistributiveInverseForProbability(T probability); }; } diff --git a/poincare/include/poincare_nodes.h b/poincare/include/poincare_nodes.h index a9b025e51..d36d674d5 100644 --- a/poincare/include/poincare_nodes.h +++ b/poincare/include/poincare_nodes.h @@ -52,6 +52,7 @@ #include #include #include +#include #include #include #include diff --git a/poincare/src/norm_cdf2.cpp b/poincare/src/norm_cdf2.cpp new file mode 100644 index 000000000..6cad822ac --- /dev/null +++ b/poincare/src/norm_cdf2.cpp @@ -0,0 +1,62 @@ +#include +#include +#include +#include +#include + +namespace Poincare { + +constexpr Expression::FunctionHelper NormCDF2::s_functionHelper; + +int NormCDF2Node::numberOfChildren() const { return NormCDF2::s_functionHelper.numberOfChildren(); } + +Expression NormCDF2Node::setSign(Sign s, ReductionContext reductionContext) { + assert(s == Sign::Positive); + return NormCDF2(this); +} + +Layout NormCDF2Node::createLayout(Preferences::PrintFloatMode floatDisplayMode, int numberOfSignificantDigits) const { + return LayoutHelper::Prefix(NormCDF2(this), floatDisplayMode, numberOfSignificantDigits, NormCDF2::s_functionHelper.name()); +} + +int NormCDF2Node::serialize(char * buffer, int bufferSize, Preferences::PrintFloatMode floatDisplayMode, int numberOfSignificantDigits) const { + return SerializationHelper::Prefix(this, buffer, bufferSize, floatDisplayMode, numberOfSignificantDigits, NormCDF2::s_functionHelper.name()); +} + +Expression NormCDF2Node::shallowReduce(ReductionContext reductionContext) { + return NormCDF2(this).shallowReduce(reductionContext); +} + +template +Evaluation NormCDF2Node::templatedApproximate(Context * context, Preferences::ComplexFormat complexFormat, Preferences::AngleUnit angleUnit) const { + Evaluation aEvaluation = childAtIndex(0)->approximate(T(), context, complexFormat, angleUnit); + Evaluation bEvaluation = childAtIndex(1)->approximate(T(), context, complexFormat, angleUnit); + Evaluation muEvaluation = childAtIndex(2)->approximate(T(), context, complexFormat, angleUnit); + Evaluation sigmaEvaluation = childAtIndex(3)->approximate(T(), context, complexFormat, angleUnit); + + T a = aEvaluation.toScalar(); + T b = bEvaluation.toScalar(); + T mu = muEvaluation.toScalar(); + T sigma = sigmaEvaluation.toScalar(); + + if (std::isnan(a) || std::isnan(b) || std::isnan(mu) || std::isnan(sigma)) { + return Complex::Undefined(); + } + if (b <= a) { + return Complex::Builder((T)0.0); + } + return Complex::Builder(NormalDistribution::CumulativeDistributiveFunctionAtAbscissa(b, mu, sigma) - NormalDistribution::CumulativeDistributiveFunctionAtAbscissa(a, mu, sigma)); +} + +Expression NormCDF2::shallowReduce(ExpressionNode::ReductionContext reductionContext) { + { + Expression e = Expression::defaultShallowReduce(); + if (e.isUndefined()) { + return e; + } + } + //TODO LEA + return *this; +} + +} diff --git a/poincare/src/normal_distribution.cpp b/poincare/src/normal_distribution.cpp index 267ef8cd4..3be5ba30e 100644 --- a/poincare/src/normal_distribution.cpp +++ b/poincare/src/normal_distribution.cpp @@ -23,8 +23,9 @@ T NormalDistribution::CumulativeDistributiveFunctionAtAbscissa(T x, T mu, T sigm return StandardNormalCumulativeDistributiveFunctionAtAbscissa((x-mu)/std::fabs(sigma)); } -double NormalDistribution::CumulativeDistributiveInverseForProbability(double probability, float mu, float sigma) { - if (sigma == 0.0f) { +template +T NormalDistribution::CumulativeDistributiveInverseForProbability(T probability, T mu, T sigma) { + if (sigma == (T)0.0) { return NAN; } return StandardNormalCumulativeDistributiveInverseForProbability(probability) * std::fabs(sigma) + mu; @@ -44,7 +45,8 @@ T NormalDistribution::StandardNormalCumulativeDistributiveFunctionAtAbscissa(T a return ((T)0.5) + ((T)0.5) * std::erf(abscissa/std::sqrt(((T)2.0))); } -double NormalDistribution::StandardNormalCumulativeDistributiveInverseForProbability(double probability) { +template +T NormalDistribution::StandardNormalCumulativeDistributiveInverseForProbability(T probability) { if (probability >= 1.0) { return INFINITY; } @@ -52,13 +54,15 @@ double NormalDistribution::StandardNormalCumulativeDistributiveInverseForProbabi return -INFINITY; } if (probability < 0.5) { - return -StandardNormalCumulativeDistributiveInverseForProbability(1-probability); + return -StandardNormalCumulativeDistributiveInverseForProbability(1.0-probability); } return std::sqrt(2.0) * erfInv(2.0 * probability - 1.0); } template float NormalDistribution::EvaluateAtAbscissa(float, float, float); -template double NormalDistribution::CumulativeDistributiveFunctionAtAbscissa(double, double, double); template float NormalDistribution::CumulativeDistributiveFunctionAtAbscissa(float, float, float); +template double NormalDistribution::CumulativeDistributiveFunctionAtAbscissa(double, double, double); +template float NormalDistribution::CumulativeDistributiveInverseForProbability(float, float, float); +template double NormalDistribution::CumulativeDistributiveInverseForProbability(double, double, double); } diff --git a/poincare/src/parsing/parser.h b/poincare/src/parsing/parser.h index 7b71cc57b..e86f738a8 100644 --- a/poincare/src/parsing/parser.h +++ b/poincare/src/parsing/parser.h @@ -121,6 +121,7 @@ private: &CommonLogarithm::s_functionHelper, &Logarithm::s_functionHelper, &NormCDF::s_functionHelper, + &NormCDF2::s_functionHelper, &PermuteCoefficient::s_functionHelper, &SimplePredictionInterval::s_functionHelper, &PredictionInterval::s_functionHelper, diff --git a/poincare/src/tree_handle.cpp b/poincare/src/tree_handle.cpp index 2f0a569a4..68e24ee6d 100644 --- a/poincare/src/tree_handle.cpp +++ b/poincare/src/tree_handle.cpp @@ -317,6 +317,7 @@ template MatrixTranspose TreeHandle::FixedArityBuilder(TreeHandle*, size_t); template NaperianLogarithm TreeHandle::FixedArityBuilder(TreeHandle*, size_t); template NormCDF TreeHandle::FixedArityBuilder(TreeHandle*, size_t); +template NormCDF2 TreeHandle::FixedArityBuilder(TreeHandle*, size_t); template NthRoot TreeHandle::FixedArityBuilder(TreeHandle*, size_t); template Opposite TreeHandle::FixedArityBuilder(TreeHandle*, size_t); template Parenthesis TreeHandle::FixedArityBuilder(TreeHandle*, size_t); diff --git a/poincare/test/approximation.cpp b/poincare/test/approximation.cpp index 993b70295..e9a64da99 100644 --- a/poincare/test/approximation.cpp +++ b/poincare/test/approximation.cpp @@ -279,6 +279,9 @@ QUIZ_CASE(poincare_approximation_function) { assert_expression_approximates_to("normcdf(1.2, 3.4, 5.6)", "0.3472125"); assert_expression_approximates_to("normcdf(1.2, 3.4, 5.6)", "3.4721249841587ᴇ-1"); + assert_expression_approximates_to("normcdf2(0.5, 3.6, 1.3, 3.4)", "0.3436388"); + assert_expression_approximates_to("normcdf2(0.5, 3.6, 1.3, 3.4)", "3.4363881299147ᴇ-1"); + assert_expression_approximates_to("permute(10, 4)", "5040"); assert_expression_approximates_to("permute(10, 4)", "5040");