From 06462490cde67826a7bc00003efa63309c18ee11 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=C3=89milie=20Feral?= Date: Mon, 15 Jan 2018 10:24:25 +0100 Subject: [PATCH] [apps] Graph: new version of the minimum search algorithm --- apps/graph/cartesian_function.cpp | 55 ++++++++++++++++--- apps/graph/cartesian_function.h | 4 +- .../graph/graph/extremum_graph_controller.cpp | 19 +++---- apps/graph/graph/extremum_graph_controller.h | 4 +- 4 files changed, 59 insertions(+), 23 deletions(-) diff --git a/apps/graph/cartesian_function.cpp b/apps/graph/cartesian_function.cpp index ebd958ca5..59f3c7810 100644 --- a/apps/graph/cartesian_function.cpp +++ b/apps/graph/cartesian_function.cpp @@ -41,13 +41,41 @@ double CartesianFunction::sumBetweenBounds(double start, double end, Poincare::C return integral.approximateToScalar(*context); } -CartesianFunction::Point CartesianFunction::mininimumBetweenBounds(double start, double end, Poincare::Context * context) const { - Point p = brentAlgorithm(start, end, context); - if (evaluateAtAbscissa(p.abscissa-k_sqrtEps, context) < p.value || evaluateAtAbscissa(p.abscissa+k_sqrtEps, context) < p.value || std::isnan(p.value)) { - p.abscissa = NAN; - p.value = NAN; +CartesianFunction::Point CartesianFunction::nextMinimumFrom(double start, double step, double max, Poincare::Context * context) const { + double bracket[3]; + Point result = {.abscissa = NAN, .value = NAN}; + double x = start; + do { + bracketMinimum(x, step, max, context, bracket); + result = brentAlgorithm(bracket[0], bracket[2], context); + x = bracket[1]; + } while (std::isnan(result.abscissa) && (step > 0.0 ? x <= max : x >= max)); + return result; +} + +void CartesianFunction::bracketMinimum(double start, double step, double max, Poincare::Context * context, double result[3]) const { + Point p[3]; + p[0] = {.abscissa = start, .value = evaluateAtAbscissa(start, context)}; + p[1] = {.abscissa = start+step, .value = evaluateAtAbscissa(start+step, context)}; + double x = start+2.0*step; + while (step > 0.0 ? x <= max : x >= max) { + p[2] = {.abscissa = x, .value = evaluateAtAbscissa(x, context)}; + if (p[0].value > p[1].value && p[2].value > p[1].value) { + result[0] = p[0].abscissa; + result[1] = p[1].abscissa; + result[2] = p[2].abscissa; + return; + } + if (p[0].value > p[1].value && p[1].value == p[2].value) { + } else { + p[0] = p[1]; + p[1] = p[2]; + } + x += step; } - return p; + result[0] = NAN; + result[1] = NAN; + result[2] = NAN; } char CartesianFunction::symbol() const { @@ -55,6 +83,9 @@ char CartesianFunction::symbol() const { } CartesianFunction::Point CartesianFunction::brentAlgorithm(double ax, double bx, Context * context) const { + if (ax > bx) { + return brentAlgorithm(bx, ax, context); + } double e = 0.0; double a = ax; double b = bx; @@ -72,9 +103,15 @@ CartesianFunction::Point CartesianFunction::brentAlgorithm(double ax, double bx, double m = 0.5*(a+b); double tol1 = k_sqrtEps*std::fabs(x)+1E-10; double tol2 = 2.0*tol1; - if (std::fabs(x-m) <= tol2-0.5*(b-a)) { - Point result = {.abscissa = x, .value = fx}; - return result; + if (std::fabs(x-m) <= tol2-0.5*(b-a)) { + double middleFax = evaluateAtAbscissa((x+a)/2.0, context); + double middleFbx = evaluateAtAbscissa((x+b)/2.0, context); + double fa = evaluateAtAbscissa(a, context); + double fb = evaluateAtAbscissa(b, context); + if (middleFax-fa <= k_sqrtEps && fx-middleFax <= k_sqrtEps && fx-middleFbx <= k_sqrtEps && middleFbx-fb <= k_sqrtEps) { + Point result = {.abscissa = x, .value = fx}; + return result; + } } double p = 0; double q = 0; diff --git a/apps/graph/cartesian_function.h b/apps/graph/cartesian_function.h index 9d5a9b91f..0593c10e9 100644 --- a/apps/graph/cartesian_function.h +++ b/apps/graph/cartesian_function.h @@ -17,12 +17,12 @@ public: double abscissa; double value; }; - Point mininimumBetweenBounds(double start, double end, Poincare::Context * context) const; + Point nextMinimumFrom(double start, double step, double max, Poincare::Context * context) const; char symbol() const override; private: constexpr static double k_sqrtEps = 1.4901161193847656E-8; // sqrt(DBL_EPSILON) constexpr static double k_goldenRatio = 0.381966011250105151795413165634361882279690820194237137864; // (3-sqrt(5))/2 - double bracketFunction(double x, double z, Poincare::Context * context); + void bracketMinimum(double start, double step, double max, Poincare::Context * context, double result[3]) const; Point brentAlgorithm(double ax, double bx, Poincare::Context * context) const; bool m_displayDerivative; }; diff --git a/apps/graph/graph/extremum_graph_controller.cpp b/apps/graph/graph/extremum_graph_controller.cpp index 869b005e2..d034decc5 100644 --- a/apps/graph/graph/extremum_graph_controller.cpp +++ b/apps/graph/graph/extremum_graph_controller.cpp @@ -24,14 +24,14 @@ View * ExtremumGraphController::view() { void ExtremumGraphController::viewWillAppear() { assert(m_function != nullptr); - CartesianFunction::Point p = computeExtremumBetweenBounds(m_graphRange->xMin(), m_graphRange->xMax()); - if (std::isnan(p.abscissa)) { + CartesianFunction::Point extremum = computeExtremumFrom(m_graphRange->xMin(), 1); + if (std::isnan(extremum.abscissa)) { m_isActive = false; m_graphView->setCursorView(nullptr); m_graphView->setBannerView(&m_defaultBannerView); } else { m_isActive = true; - m_cursor->moveTo(p.abscissa, p.value); + m_cursor->moveTo(extremum.abscissa, extremum.value); m_graphRange->panToMakePointVisible(m_cursor->x(), m_cursor->y(), k_cursorTopMarginRatio, SimpleInteractiveCurveViewController::k_cursorRightMarginRatio, k_cursorBottomMarginRatio, SimpleInteractiveCurveViewController::k_cursorLeftMarginRatio); reloadBannerView(); } @@ -67,11 +67,7 @@ void ExtremumGraphController::reloadBannerView() { } bool ExtremumGraphController::moveCursor(int direction) { - double x = m_cursor->x(); - float step = m_graphRange->xGridUnit()/SimpleInteractiveCurveViewController::k_numberOfCursorStepsInGradUnit; - float start = direction > 0 ? x+step : m_graphRange->xMin(); - float end = direction > 0 ? m_graphRange->xMax() : x-step; - CartesianFunction::Point newExtremum = computeExtremumBetweenBounds(start, end); + CartesianFunction::Point newExtremum = computeExtremumFrom(m_cursor->x(), direction); if (std::isnan(newExtremum.abscissa)) { return false; } @@ -89,9 +85,12 @@ const char * MinimumGraphController::title() { return I18n::translate(I18n::Message::Minimum); } -CartesianFunction::Point MinimumGraphController::computeExtremumBetweenBounds(float start, float end) { +CartesianFunction::Point MinimumGraphController::computeExtremumFrom(double start, int direction) { TextFieldDelegateApp * myApp = (TextFieldDelegateApp *)app(); - return m_function->mininimumBetweenBounds(start, end, myApp->localContext()); + double step = m_graphRange->xGridUnit()/100.0; + step = direction < 0 ? -step : step; + double max = direction > 0 ? m_graphRange->xMax() : m_graphRange->xMin(); + return m_function->nextMinimumFrom(start, step, max, myApp->localContext()); } } diff --git a/apps/graph/graph/extremum_graph_controller.h b/apps/graph/graph/extremum_graph_controller.h index aa07e8f21..72c932ea0 100644 --- a/apps/graph/graph/extremum_graph_controller.h +++ b/apps/graph/graph/extremum_graph_controller.h @@ -23,7 +23,7 @@ protected: BannerView * bannerView() override { return m_bannerView; } void reloadBannerView(); bool moveCursor(int direction); - virtual CartesianFunction::Point computeExtremumBetweenBounds(float start, float end) = 0; + virtual CartesianFunction::Point computeExtremumFrom(double start, int direction) = 0; GraphView * m_graphView; BannerView * m_bannerView; Shared::InteractiveCurveViewRange * m_graphRange; @@ -38,7 +38,7 @@ public: MinimumGraphController(Responder * parentResponder, GraphView * graphView, BannerView * bannerView, Shared::InteractiveCurveViewRange * curveViewRange, Shared::CurveViewCursor * cursor); const char * title() override; private: - CartesianFunction::Point computeExtremumBetweenBounds(float start, float end) override; + CartesianFunction::Point computeExtremumFrom(double start, int direction) override; }; }