Every step a track takes is at most max_step_size, the first one included

CHANGELOG.md

@@ -77,6 +77,7 @@ Bug fixes for 4.0.
 - Differentiating an expression differentiates each distinct subexpression once, however many paths in the expression reach it (#505).  It used to differentiate a shared subexpression once per path, so a determinant built by cofactor expansion with its sub-minors shared took time proportional to k! -- 24 s for the three partial derivatives of a 7 x 7 determinant of quadratics, all of it paid on a `System`'s first evaluation, which builds its Jacobian.  `Node::Differentiate` now keeps each derivative for the length of the outermost call, and a Jacobian is built under one `node::DerivativeScope`, so subexpressions its functions share are differentiated once among all of them.  A node type's own rule is now `DifferentiateImpl`; a Python subclass of `Node` still overrides `Differentiate`.
 - Building a sum or product no longer walks its operands (#505).  Putting the operands in canonical order asks each one's variables, multidegree and, on a tie, canonical encoding, and each answer was a walk of the operand's DAG -- of its whole expansion, for `Degree(v)` under a negation -- so building an expression cost the square of its size, or worse.  Each node now keeps a degree record, made once from its operands' records: its variables, its degree and multidegree in each variable it mentions, its total degree, and its answer for every variable it does not mention.  `Degree`, `MultiDegree` and `GatherVariables` read it.  On the Barth sextic's 50 deflated critical-point systems, the first evaluations took 121 s and take 0.6 s.  A node type's degree rule is now `DegreeImpl`; a Python subclass of `Node` still overrides `Degree`.
 - **The canonicalization tie-break is a content key, and the system encoding version is `b2sysenc/3`** (#505, ADR-0073).  Operands that tie on multidegree were ordered by their whole canonical encoding, which costs a walk of the operand.  They are now ordered by a content key: the SHA-256 of the operand's own canonical encoding with each of its operands written as that operand's key -- a Merkle hash, made from the operands' keys once per node.  The first evaluation of a 7 x 7 determinant of quadratics took 24.6 s in 4.0.0 and takes 0.04 s.  Where operands tie, their order in a system's canonical encoding changes, so the encoding version moves to `b2sysenc/3` and every system digest changes: records written by 4.0.0 are asked again rather than recalled, and `System.from_canonical` refuses `b2sysenc/2` text.
+- Every step a track takes is at most `max_step_size`, the first one included (#508).  A track opened with `initial_step_size` whatever `max_step_size` said, and `max_step_size` bounded only the step's growth after that, so with `max_step_size` 1e-3 and the default `initial_step_size` 1/10 the adaptive tracker took 17 steps above the maximum before reaching it: six of 0.1, six of 0.05, then halving down.  A track asked to cross a stretch of the path in small steps could cross it in one, and land on another path.  A track now opens with the smaller of the two (`Tracker::FirstStepSize`).  With the defaults, where both are 1/10, nothing changes.
 
 _______________________________________________________________________________
 

core/include/bertini2/trackers/amp_tracker.hpp

@@ -549,7 +549,7 @@ namespace bertini{
                 if (reinitialize_stepsize_)
                 {
                     real_mp segment_length = abs(start_time-end_time)/Get<Stepping>().min_num_steps;
-                    SetStepSize(min(real_mp(Get<Stepping>().initial_step_size, initial_precision_),segment_length));
+                    SetStepSize(min(real_mp(FirstStepSize(), initial_precision_),segment_length));
                 }
                 current_stepsize_.precision(initial_precision_);
 

core/include/bertini2/trackers/base_tracker.hpp

@@ -206,7 +206,7 @@ namespace bertini{
                 this->template Set<SteppingConfig>(stepping);
                 this->template Set<NewtonConfig>(newton);
 
-                current_stepsize_ = BaseRealT(stepping.initial_step_size);
+                current_stepsize_ = BaseRealT(FirstStepSize());
             }
 
 
@@ -598,6 +598,20 @@ namespace bertini{
             }
 
 
+            /**
+            \brief The length of a track's first step, before any limit from the length of the track.
+
+            The stepping settings' initial_step_size, but never more than their max_step_size: every
+            step a track takes is at most max_step_size, the first one included (#508).
+            */
+            auto FirstStepSize() const
+            {
+                auto const& stepping = this->template Get<Stepping>();
+                return stepping.max_step_size < stepping.initial_step_size ? stepping.max_step_size
+                                                                           : stepping.initial_step_size;
+            }
+
+
             /**
             \brief Switch resetting of initial step size to that of the stepping settings.
 

core/include/bertini2/trackers/fixed_precision_tracker.hpp

@@ -486,7 +486,7 @@ namespace bertini{
                 this->endtime_ = end_time;
                 std::get<Vec<BaseComplexT> >(this->current_space_) = start_point;
                 if (this->reinitialize_stepsize_)
-                    this->SetStepSize(min(BaseRealT(Get<Stepping>().initial_step_size),abs(start_time-end_time)/Get<Stepping>().min_num_steps));
+                    this->SetStepSize(min(BaseRealT(this->FirstStepSize()),abs(start_time-end_time)/Get<Stepping>().min_num_steps));
 
                 ResetCounters();
 
@@ -615,7 +615,7 @@ namespace bertini{
                 this->endtime_ = end_time;
                 std::get<Vec<BaseComplexT> >(this->current_space_) = start_point;
                 if (this->reinitialize_stepsize_)
-                    this->SetStepSize(min(real_mp(Get<Stepping>().initial_step_size),real_mp(abs(start_time-end_time)/Get<Stepping>().min_num_steps)));
+                    this->SetStepSize(min(real_mp(this->FirstStepSize()),real_mp(abs(start_time-end_time)/Get<Stepping>().min_num_steps)));
 
                 ResetCounters();
 

core/test/tracking_basics/amp_tracker_test.cpp

@@ -1753,4 +1753,67 @@ BOOST_AUTO_TEST_CASE(a_tracker_never_returns_never_started)
     BOOST_CHECK(distinct.size() > 1);                                // and at least one did not
 }
 
+
+namespace {
+// The time after each successful step, as a double.
+struct StepTimeRecorder : public bertini::Observer<bertini::tracking::AMPTracker>
+{
+    using Emitter = bertini::tracking::TrackerTraits<bertini::tracking::AMPTracker>::EventEmitterType;
+
+    std::vector<double> times; ///< The real part of the time after each successful step.
+
+    bertini::ObserveResult Observe(bertini::AnyEvent const& e) override
+    {
+        if (auto p = dynamic_cast<const bertini::tracking::SuccessfulStep<Emitter>*>(&e))
+            times.push_back(static_cast<double>(real(p->Get().CurrentTime())));
+        return bertini::ObserveResult::KeepObserving;
+    }
+};
+}
+
+
+/**
+Every step a track takes is at most max_step_size, the first one included (b2#508).
+
+A track opened with initial_step_size whatever max_step_size said, and max_step_size bounded only
+the step's growth after that: a track asked for steps of at most 1/1000 opened with the default
+1/10, crossing in one step the stretch its caller had asked it to cross in a hundred.
+*/
+BOOST_AUTO_TEST_CASE(every_step_is_at_most_the_max_step_size_the_first_included)
+{
+    using namespace bertini::tracking;
+    using Cmp = bertini::complex_mp;
+    DefaultPrecision(30);
+
+    auto y = Variable::Make("y");
+    auto t = Variable::Make("t");
+    System sys;
+    sys.AddFunction(y - t);
+    sys.AddPathVariable(t);
+    sys.AddVariableGroup(VariableGroup{y});
+
+    AMPTracker tracker(sys);
+    SteppingConfig stepping;
+    stepping.max_step_size = bertini::mpq_rational(1, 1000);
+    stepping.max_num_steps = 10000;
+    tracker.Setup(Predictor::RKF45, 1e-5, 1e5, stepping, NewtonConfig());
+    tracker.PrecisionSetup(AMPConfigFrom(sys));
+
+    StepTimeRecorder recorder;
+    tracker.AddObserver(recorder);
+    bertini::Vec<Cmp> start(1), end_point;
+    start << Cmp(1);
+    auto code = tracker.TrackPath(end_point, Cmp(1), Cmp(0), start);
+    tracker.RemoveObserver(recorder);
+
+    BOOST_CHECK(code == bertini::SuccessCode::Success);
+    BOOST_REQUIRE(!recorder.times.empty());
+    double before = 1;
+    for (double after : recorder.times)
+    {
+        BOOST_CHECK_LE(before - after, 1e-3 * (1 + 1e-12));
+        before = after;
+    }
+}
+
 BOOST_AUTO_TEST_SUITE_END()

core/test/tracking_basics/fixed_precision_tracker_test.cpp

@@ -523,4 +523,69 @@ BOOST_AUTO_TEST_CASE(a_truncation_threshold_set_through_the_config_truncates)
 }
 
 
+namespace {
+// The time after each successful step, as a double.
+template<typename TrackerT>
+struct StepTimeRecorder : public bertini::Observer<TrackerT>
+{
+    using Emitter = typename bertini::tracking::TrackerTraits<TrackerT>::EventEmitterType;
+
+    std::vector<double> times; ///< The real part of the time after each successful step.
+
+    bertini::ObserveResult Observe(bertini::AnyEvent const& e) override
+    {
+        if (auto p = dynamic_cast<const bertini::tracking::SuccessfulStep<Emitter>*>(&e))
+            times.push_back(static_cast<double>(real(p->Get().CurrentTime())));
+        return bertini::ObserveResult::KeepObserving;
+    }
+};
+
+// Tracks y = t from t = 1 to 0 with max_step_size 1/1000 and the default initial_step_size,
+// 1/10, and returns the time after each step.
+template<typename TrackerT, typename ComplexT>
+std::vector<double> StepTimesUnderAMaxStepOfOneThousandth()
+{
+    using namespace bertini::tracking;
+
+    Var y = Variable::Make("y");
+    Var t = Variable::Make("t");
+    System sys;
+    sys.AddFunction(y - t);
+    sys.AddPathVariable(t);
+    sys.AddVariableGroup(VariableGroup{y});
+
+    TrackerT tracker(sys);
+    SteppingConfig stepping;
+    stepping.max_step_size = bertini::mpq_rational(1, 1000);
+    stepping.max_num_steps = 10000;
+    tracker.Setup(Predictor::RKF45, 1e-5, 1e5, stepping, NewtonConfig());
+
+    StepTimeRecorder<TrackerT> recorder;
+    tracker.AddObserver(recorder);
+    Vec<ComplexT> start(1), end_point;
+    start << ComplexT(1);
+    tracker.TrackPath(end_point, ComplexT(1), ComplexT(0), start);
+    tracker.RemoveObserver(recorder);
+    return recorder.times;
+}
+}
+
+
+// b2#508: a track's first step was initial_step_size whatever max_step_size said, so a track
+// asked to take steps of at most 1/1000 opened with a step of 1/10.
+BOOST_AUTO_TEST_CASE(the_first_step_is_at_most_the_max_step_size)
+{
+    DefaultPrecision(30);
+    using namespace bertini::tracking;
+
+    auto const in_double = StepTimesUnderAMaxStepOfOneThousandth<DoublePrecisionTracker, complex_dbl>();
+    auto const in_multiple = StepTimesUnderAMaxStepOfOneThousandth<MultiplePrecisionTracker, mpfr>();
+    for (auto const& times : {in_double, in_multiple})
+    {
+        BOOST_REQUIRE(!times.empty());
+        BOOST_CHECK_LE(1 - times.front(), 1e-3 * (1 + 1e-12));
+    }
+}
+
+
 BOOST_AUTO_TEST_SUITE_END()