cauchy.hpp

Go to the documentation of this file.

//This file is part of Bertini 2.
//
//cauchy_endgame.hpp is free software: you can redistribute it and/or modify
//it under the terms of the GNU General Public License as published by
//the Free Software Foundation, either version 3 of the License, or
//(at your option) any later version.
//
//cauchy_endgame.hpp is distributed in the hope that it will be useful,
//but WITHOUT ANY WARRANTY; without even the implied warranty of
//MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//GNU General Public License for more details.
//
//You should have received a copy of the GNU General Public License
//along with cauchy_endgame.hpp.  If not, see <http://www.gnu.org/licenses/>.
//
// Copyright(C) Bertini2 Development Team
//
// See <http://www.gnu.org/licenses/> for a copy of the license,
// as well as COPYING.  Bertini2 is provided with permitted
// additional terms in the b2/licenses/ directory.
 
#pragma once
 
#include "bertini2/endgames/base_endgame.hpp"
#include "bertini2/trackers/adaptive_precision_utilities.hpp"  // tracking::adaptive::SetPrecision, for container migration
 

 
namespace bertini{ namespace endgame{
 

template<typename PrecT>
class CauchyEndgame :
    public virtual EndgameBase<CauchyEndgame<PrecT>, PrecT>
{
public:
    using BaseEGT = EndgameBase<CauchyEndgame<PrecT>, PrecT>;  
    using FinalEGT = CauchyEndgame<PrecT>;  
    using TrackerType = typename PrecT::TrackerType;  
 
    using BaseComplexT = typename tracking::TrackerTraits<TrackerType>::BaseComplexT;  
    using BaseRealT = typename tracking::TrackerTraits<TrackerType>::BaseRealT;  
 
    using EmitterType = CauchyEndgame<PrecT>;  
 
protected:
 
    using EndgameBase<CauchyEndgame<PrecT>, PrecT>::NotifyObservers;
 
 
 
 
    using TupleOfTimes = typename BaseEGT::TupleOfTimes;  
    using TupleOfSamps = typename BaseEGT::TupleOfSamps;  
    using TupOfVec = typename BaseEGT::TupOfVec;  
 
    using BCT = BaseComplexT;  
    using BRT = BaseRealT;  
 
    using Configs = typename AlgoTraits<FinalEGT>::NeededConfigs;  
    using ConfigsAsTuple = typename Configs::ToTuple;  

    mutable TupleOfTimes pseg_times_;
    mutable TupleOfSamps pseg_samples_; //samples used for the first approximation.
    mutable TupleOfTimes cauchy_times_;
    mutable TupleOfSamps cauchy_samples_;
 
    // The c/k probe lives on the base endgame now: the operating-zone test is shared with the
    // power series endgame (b2#402).
 
    // Scratch for LatestTimeImpl to return a BCT reference when the endgame is computing in the
    // complex_dbl fast lane (the latest time then lives in the complex_dbl slot, not the BCT slot).
    mutable BCT latest_time_cache_; 
 
 
 
 
 
public:
 
 
 

    template<typename ComplexT>
    void ClearTimesAndSamples()
    {
        std::get<TimeCont<ComplexT> >(pseg_times_).clear();
        std::get<TimeCont<ComplexT> >(cauchy_times_).clear();
        std::get<SampCont<ComplexT> >(pseg_samples_).clear();
        std::get<SampCont<ComplexT> >(cauchy_samples_).clear();}

    template<typename ComplexT>
    void SetPSEGTimes(TimeCont<ComplexT> pseg_times_to_set)
    { std::get<TimeCont<ComplexT> >(pseg_times_) = pseg_times_to_set;}

    template<typename ComplexT>
    TimeCont<ComplexT>& GetPSEGTimes() {return std::get<TimeCont<ComplexT> >(pseg_times_);}
    template<typename ComplexT>
    const TimeCont<ComplexT>& GetPSEGTimes() const {return std::get<TimeCont<ComplexT> >(pseg_times_);}

    template<typename ComplexT>
    void SetPSEGSamples(SampCont<ComplexT> const& pseg_samples_to_set) { std::get<SampCont<ComplexT> >(pseg_samples_) = pseg_samples_to_set;}

    template<typename ComplexT>
    SampCont<ComplexT>& GetPSEGSamples() {return std::get<SampCont<ComplexT> >(pseg_samples_);}
    template<typename ComplexT>
    const SampCont<ComplexT>& GetPSEGSamples() const {return std::get<SampCont<ComplexT> >(pseg_samples_);}
    template<typename ComplexT>
    void SetCauchySamples(SampCont<ComplexT> const& cauchy_samples_to_set)
    {
        std::get<SampCont<ComplexT> >(cauchy_samples_) = cauchy_samples_to_set;
    }

    template<typename ComplexT>
    SampCont<ComplexT>& GetCauchySamples()
    {
        return std::get<SampCont<ComplexT> >(cauchy_samples_);
    }

    template<typename ComplexT>
    const SampCont<ComplexT>& GetCauchySamples() const { return std::get<SampCont<ComplexT> >(cauchy_samples_); }
 

    template<typename ComplexT>
    void SetCauchyTimes(TimeCont<ComplexT> const& cauchy_times_to_set)
    {
        std::get<TimeCont<ComplexT> >(cauchy_times_) = cauchy_times_to_set;
    }

    template<typename ComplexT>
    TimeCont<ComplexT>& GetCauchyTimes()
    {
        return std::get<TimeCont<ComplexT> >(cauchy_times_);
    }

    template<typename ComplexT>
    const TimeCont<ComplexT>& GetCauchyTimes() const
    {
        return std::get<TimeCont<ComplexT> >(cauchy_times_);
    }
 

    const BCT& LatestTimeImpl() const
    {
        // In the adaptive-numeric-type endgame the latest time may live in the complex_dbl slot (the
        // fast lane), with the BCT slot empty.  Dispatch on which slot actually holds data, so this is
        // correct both during the run (observer events) and after it (solution metadata), regardless of
        // the adaptive_numeric_type_active_ flag.  Fixed precision compiles to the original BCT read.
        if constexpr (tracking::TrackerTraits<TrackerType>::IsAdaptivePrec)
        {
            if (GetPSEGTimes<BCT>().empty())
            {
                auto const& dbl_times = GetPSEGTimes<complex_dbl>();
                latest_time_cache_ = dbl_times.empty() ? BCT(0) : BCT(dbl_times.back());
                return latest_time_cache_;
            }
        }
        return GetPSEGTimes<BCT>().back();
    }
 

    void SetCauchySettings(CauchyConfig const& new_cauchy_settings)
    {
        this->template Set<CauchyConfig>(new_cauchy_settings);
    }

    const auto& GetCauchySettings() const
    {
        return this->template Get<CauchyConfig>();
    }
 

    explicit CauchyEndgame(TrackerType const& tr,
                            const ConfigsAsTuple& settings )
      : EndgamePrecPolicyBase<TrackerType>(tr), BaseEGT(tr, settings)
    { }

    template< typename... Ts >
        CauchyEndgame(TrackerType const& tr, const Ts&... ts ) : CauchyEndgame(tr, Configs::Unpermute( ts... ) )
        {}
 
 
    virtual ~CauchyEndgame() = default;

    void ValidateConfigs()
    {
        if (this->EndgameSettings().num_sample_points < 3) // need to make sure we won't track right through the origin.
        {
            std::stringstream err_msg;
            err_msg << "ERROR: The number of sample points " << this->EndgameSettings().num_sample_points << " for circle tracking must be >= 3";
            throw std::runtime_error(err_msg.str());
        }
    }

    template<typename ComplexT>
    SuccessCode CircleTrack(ComplexT const& target_time)
    {
        using bertini::Precision;
        using RealT = typename Eigen::NumTraits<ComplexT>::Real;
        using std::acos;
 
        ValidateConfigs();
 
        auto& circle_times = std::get<TimeCont<ComplexT> >(cauchy_times_);
        auto& circle_samples = std::get<SampCont<ComplexT> >(cauchy_samples_);
 
        ComplexT starting_time = circle_times.back();  // take a COPY here, so won't invalidate it later
 
        // the initial sample has already been added to the sample repo... so don't do that here, please
 
        const auto num_vars = this->GetSystem().NumVariables();
 
        for (unsigned ii = 0; ii < this->EndgameSettings().num_sample_points; ++ii)
        {
            const Vec<ComplexT>& current_sample = circle_samples.back();
            const ComplexT& current_time = circle_times.back();
 
#ifndef BERTINI_DISABLE_PRECISION_CHECKS
            if (Precision(current_time)!=Precision(current_sample)){
                std::stringstream err_msg;
                err_msg << "current time and sample for circle track must be of same precision.  respective precisions: " << Precision(current_time) << " " << Precision(current_sample) << std::endl;
                throw std::runtime_error(err_msg.str());
            }
#endif
 
            //set up the time value for the next sample.
            using std::polar;
 
#ifndef USE_BMP_COMPLEX
            using bertini::polar;
#endif
 
            //Generalized since we could have a nonzero target time.
            using std::arg;
            RealT radius = abs(starting_time - target_time), angle = arg(starting_time - target_time); // generalized for nonzero target_time.
 
            auto next_sample = Vec<ComplexT>(num_vars);
            ComplexT next_time = (ii==this->EndgameSettings().num_sample_points-1)
                                ?
                              starting_time
                                :
                              polar(radius, (ii+1)*2*acos(static_cast<RealT>(-1)) / (this->EndgameSettings().num_sample_points) + angle) + target_time;
            // If we are tracking to a nonzero target time we need to shift our values to track to. This is a step that may not be needed if target_time = 0
                              ;
 
 
            auto tracking_success = this->EndgameTrackPath(next_sample, current_time, next_time, current_sample);
            if (tracking_success != SuccessCode::Success)
            {
                return tracking_success;
            }
 
            // Pure-(i) numeric-type escalation: if the tracker's authoritative precision climbed above
            // the precision this endgame is computing in, double no longer suffices for this circle.
            // Bail to the migrate-and-retry driver (RunImplAMP), which discards this partial circle,
            // migrates the durable state to mpfr, and re-tracks the circle in mpfr.  Compile-time elided
            // for fixed precision (and shields GetCurrentPrecision(), which fixed trackers lack).
            if constexpr (tracking::TrackerTraits<TrackerType>::IsAdaptivePrec)
                if (this->adaptive_numeric_type_active_ &&
                    this->GetTracker().GetCurrentPrecision() > this->current_endgame_precision_)
                    return SuccessCode::HigherPrecisionNecessary;
 
            // CircleAdvanced carries the new point/time at BaseComplexT.  The fixed/mpfr lane emits exactly
            // as before.  The complex_dbl fast lane would have to convert the point to mpfr for the event,
            // so we only pay that when something is actually observing (temporaries live through the
            // synchronous NotifyObservers).
            if constexpr (std::is_same<ComplexT, BCT>::value)
                NotifyObservers(CircleAdvanced<EmitterType>(*this, next_sample, next_time));
            else if (this->HasObservers())
            {
                Vec<BCT> ev_pt(next_sample.size());
                for (Eigen::Index i = 0; i < next_sample.size(); ++i) ev_pt(i) = BCT(next_sample(i));
                BCT ev_t(next_time);
                NotifyObservers(CircleAdvanced<EmitterType>(*this, ev_pt, ev_t));
            }
 
 
            this->EnsureAtPrecision(next_time,Precision(next_sample)); assert(Precision(next_time)==Precision(next_sample));
 
            // auto refinement_success = this->RefineSample(next_sample, next_sample, next_time,
            //                          this->FinalTolerance() * this->EndgameSettings().sample_point_refinement_factor,
            //                          this->EndgameSettings().max_num_refinements);
            // if (refinement_success != SuccessCode::Success)
            // {
            //  return refinement_success;
            // }
 
            // this->EnsureAtPrecision(next_time,Precision(next_sample));
            // assert(Precision(next_time)==Precision(next_sample));
 
            AddToCauchyData(next_time, next_sample);
 
            // NotifyObservers(SampleRefined<EmitterType>(*this));
        }
 
        return SuccessCode::Success;
 
    }//end CircleTrack

    template<typename ComplexT>
    void AddToCauchyData(ComplexT const& time, Vec<ComplexT> const& sample)
    {
        std::get<TimeCont<ComplexT>>(cauchy_times_).push_back(time);
        std::get<SampCont<ComplexT>>(cauchy_samples_).push_back(sample);
    }

    template<typename ComplexT>
    void AddToPSData(ComplexT const& time, Vec<ComplexT> const& sample)
    {
        std::get<TimeCont<ComplexT>>(pseg_times_).push_back(time);
        std::get<SampCont<ComplexT>>(pseg_samples_).push_back(sample);
    }
 
    // GetCOverKProbe, ComputeCOverK and CheckForCOverKStabilization now live on the base endgame:
    // the estimate reads only the geometric approach samples, which every flavor maintains, so the
    // operating-zone test is shared rather than Cauchy's alone (b2#402).  The Cauchy endgame
    // passes its own pseg_samples_ to ComputeCOverK.
 
 
 
    /*
    Input: A time value and the space value above that time.
 
    Output: An real_mp representing a tolerance threshold for declaring a loop to be closed.
    Details: Used in Bertini 1 as a heuristic for computing separatedness of roots. Decided to not be used since assumptions for this tolerance are not usually met.
    template<typename ComplexT>
    real_mp FindToleranceForClosedLoop(ComplexT x_time, Vec<ComplexT> x_sample)
    {
        auto degree_max = std::max(this->GetTracker().AMP_config_.degree_bound,real_mp("2.0"));
        auto K = this->GetTracker().AMP_config_.coefficient_bound;
        real_mp N;
        real_mp M;
        real_mp L;
        if(max_closed_loop_tolerance_ < min_closed_loop_tolerance_)
        {
            max_closed_loop_tolerance_ = min_closed_loop_tolerance_;
        }
        auto error_tolerance = real_mp("1e-13");
        if(x_sample.size() <= 1)
        {
            N = degree_max;
        }
        else
        {
            N = ComputeCombination(degree_max + x_sample[0].precision() - 1, x_sample[0].precision() - 1);
        }
        M = degree_max * (degree_max - 1) * N;
        auto jacobian_at_current_time = this->GetSystem().Jacobian(x_sample,x_time);
        auto minimum_singular_value = Eigen::JacobiSVD< Mat<ComplexT> >(jacobian_at_current_time).singularValues()(this->GetSystem().NumVariables() - 1 );
        auto norm_of_sample = x_sample.norm();
        L = pow(norm_of_sample,degree_max - 2);
        auto tol = K * L * M;
        if (tol == 0) // fail-safe
                    tol = minimum_singular_value;
            else
            {
            tol = real_mp("2.0") / tol;
        tol = tol * minimum_singular_value;
            }
            // make sure that tol is between min_tol & max_tol
        if (tol > max_closed_loop_tolerance_) // tol needs to be <= max_tol
            tol = max_closed_loop_tolerance_;
        if (tol < min_closed_loop_tolerance_) // tol needs to be >= min_tol
            tol = min_closed_loop_tolerance_;
            return tol;
    }// end FindToleranceForClosedLoop
    */
 
 

    template<typename ComplexT>
    bool CheckClosedLoop()
    {
        auto& times = std::get<TimeCont<ComplexT> >(cauchy_times_);
        auto& samples = std::get<SampCont<ComplexT> >(cauchy_samples_);
 
        if((samples.front() - samples.back()).template lpNorm<Eigen::Infinity>() < this->GetTracker().TrackingTolerance())
        {
            return true;
        }
 
        if constexpr (tracking::TrackerTraits<TrackerType>::IsAdaptivePrec)
        {
            //Ensure all samples are of the same precision.
            auto new_precision = this->EnsureAtUniformPrecision(times, samples);
        }
 
        this->GetTracker().Refine(samples.front(),samples.front(),times.front(),this->FinalTolerance(),this->EndgameSettings().max_num_refinements);
        this->GetTracker().Refine(samples.back(),samples.back(),times.back(),this->FinalTolerance(),this->EndgameSettings().max_num_refinements);
 
        if((samples.front() - samples.back()).template lpNorm<Eigen::Infinity>() < this->GetTracker().TrackingTolerance())
        {
            return true;
        }
        return false;
 
    }//end CheckClosedLoop
 
 
 

    template<typename ComplexT>
    bool RatioEGOperatingZoneTest(ComplexT const& target_time) const
    {
        using RealT = typename Eigen::NumTraits<ComplexT>::Real;
        RealT min(1e300);
        RealT max(0);
        auto& times = std::get<TimeCont<ComplexT> >(cauchy_times_);
        auto& samples = std::get<SampCont<ComplexT> >(cauchy_samples_);
        if(norm(times.front() - target_time) < GetCauchySettings().ratio_cutoff_time)
        {
            return true;
        }
        else
        {
            RealT norm;
            for(unsigned int ii=0; ii < this->EndgameSettings().num_sample_points; ++ii)
            {
                norm = samples[ii].template lpNorm<Eigen::Infinity>();
                if(norm > max)
                {
                    max = norm;
                }
                if(norm < min)
                {
                    min = norm;
                }
            }
 
            if(min > this->FinalTolerance() && max > this->FinalTolerance())
            {
                norm = min / max;
                if(norm < GetCauchySettings().maximum_cauchy_ratio && (max - min) > this->FinalTolerance())
                {
                    return false; // bad ratio and too far apart.
                }
            }
        }
        return true;
 
    }//end RatioEGOperatingZoneTest
 
 
    /*
    Input: All input needed is available as class data members.
    Output: SuccessCode declaring if we are successful or not.
    Details: Computes cauchy samples using Circle Track. Compares the ratios of the maximum and minimum norms of cauchy samples
             using RatioEGOperatingZoneTest. Then attempts to see if we have closed our loop, otherwise it continues the process outlined.
    */
 
 

    template<typename ComplexT>
    SuccessCode InitialCauchyLoops(ComplexT const& target_time)
    {
        using std::max;
        // auto fail_safe_max_cycle_number = max(GetCauchySettings().fail_safe_maximum_cycle_number,this->CycleNumber());
        auto fail_safe_max_cycle_number = GetCauchySettings().fail_safe_maximum_cycle_number;
 
        auto initial_cauchy_loop_success = SuccessCode::Success;
 
        bool loop_hasnt_closed = true;
        while (loop_hasnt_closed)
        {
            this->CycleNumber(0);
            ClearAndSeedCauchyData<ComplexT>();
 
            // track around a circle once.  we'll use it to measure whether we believe we are in the eg operating zone, based on the ratio of norms of sample points around the circle
            auto tracking_success = CircleTrack(target_time);
            this->IncrementCycleNumber(1);
            if (tracking_success != SuccessCode::Success)
                return tracking_success;
 
            // find the ratio of the maximum and minimum coordinate wise for the loop.
            if (RatioEGOperatingZoneTest<ComplexT>(target_time))
            { // then we believe we are in the EG operating zone, since the path is relatively flat.  i still disbelieve this is a good test (dab 20160310)
                while (true)
                {
                    if (CheckClosedLoop<ComplexT>())
                    {//error is small enough, exit the loop with success.
                        NotifyObservers(ClosedLoop<EmitterType>(*this));
                        initial_cauchy_loop_success = SuccessCode::Success;
                        loop_hasnt_closed = false;
                        break;
                    }
                    else if(this->CycleNumber() > fail_safe_max_cycle_number)
                    {// too many iterations
                        initial_cauchy_loop_success = SuccessCode::CycleNumTooHigh;
                        loop_hasnt_closed = false;
                        break;
                    }
 
                    //compute next loop, the last sample in times and samples is the sample our loop ended on. Either where we started or on another sheet at the same time value.
                    tracking_success = CircleTrack(target_time);
                    this->IncrementCycleNumber(1);
                    if(tracking_success != SuccessCode::Success)
                        return tracking_success;
                }
            }//end if (RatioEGOperatingZoneTest())
            else
            {
                auto advance_success = AdvanceTime<ComplexT>(target_time);
                if (advance_success!=SuccessCode::Success)
                    return advance_success;
            }
        } //end while(loop_hasnt_closed)
 
        return initial_cauchy_loop_success;
    }//end InitialCauchyLoops
 
 

    template <typename ComplexT>
    void RotateOntoPS(ComplexT const& next_time, Vec<ComplexT> const& next_sample)
    {
        auto& ps_times = std::get<TimeCont<ComplexT> >(pseg_times_);
        auto& ps_samples = std::get<SampCont<ComplexT> >(pseg_samples_);
 
        ps_times.pop_front();
        ps_samples.pop_front();
 
        ps_times.push_back(next_time);
        ps_samples.push_back(next_sample);
    }

    template <typename ComplexT>
    void ClearAndSeedCauchyData()
    {
        auto& cau_times = std::get<TimeCont<ComplexT> >(cauchy_times_);
        auto& cau_samples = std::get<SampCont<ComplexT> >(cauchy_samples_);
        auto& ps_times = std::get<TimeCont<ComplexT> >(pseg_times_);
        auto& ps_samples = std::get<SampCont<ComplexT> >(pseg_samples_);
 
        cau_times.clear();
        cau_samples.clear();
        cau_samples.push_back(ps_samples.back());
        cau_times.push_back(ps_times.back());
    }
 

    template<typename ComplexT>
    SuccessCode InitialApproximation(ComplexT const& start_time, Vec<ComplexT> const& start_point,
                                                ComplexT const& target_time, Vec<ComplexT> & approximation)
    {
        auto init_success = GetIntoEGZone(start_time, start_point, target_time);
        if (init_success!= SuccessCode::Success)
            return init_success;
 
        auto cauchy_loop_success = InitialCauchyLoops<ComplexT>(target_time);
        if (cauchy_loop_success != SuccessCode::Success)
            return cauchy_loop_success;
 
        return ComputeCauchyApproximationOfXAtT0(approximation);
 
    }//end InitialApproximation
 
 

    template<typename ComplexT>
    SuccessCode GetIntoEGZone(ComplexT const& start_time, Vec<ComplexT> const& start_point, ComplexT const& target_time)
    {
        using RealT = typename Eigen::NumTraits<ComplexT>::Real;
 
        //initialize array holding c_over_k estimates
        std::deque<RealT> c_over_k;
 
        auto& ps_times = std::get<TimeCont<ComplexT> >(pseg_times_);
        auto& ps_samples = std::get<SampCont<ComplexT> >(pseg_samples_);
 
        //Compute initial samples for pseg
        auto initial_sample_success = this->ComputeInitialSamples(start_time, target_time, start_point, ps_times, ps_samples);
        if (initial_sample_success!=SuccessCode::Success)
            return initial_sample_success;
 
        c_over_k.push_back(this->template ComputeCOverK<ComplexT>(std::get<SampCont<ComplexT>>(pseg_samples_)));
 
 
        //track until for more c_over_k estimates or until we reach a cutoff time.
        for (unsigned ii = 0; ii < this->EndgameSettings().num_needed_for_stabilization; ++ii)
        {
            auto advance_success = AdvanceTime<ComplexT>(target_time);
            if (advance_success!=SuccessCode::Success)
                return advance_success;
            c_over_k.push_back(this->template ComputeCOverK<ComplexT>(std::get<SampCont<ComplexT>>(pseg_samples_)));
        }//end while
 
 
        //have we stabilized yet?
        while(!this->CheckForCOverKStabilization(c_over_k) && abs(ps_times.back()-target_time) > GetCauchySettings().cycle_cutoff_time)
        {
            auto advance_success = AdvanceTime<ComplexT>(target_time);
            if (advance_success!=SuccessCode::Success)
                return advance_success;
 
            c_over_k.pop_front();
            c_over_k.push_back(this->template ComputeCOverK<ComplexT>(std::get<SampCont<ComplexT>>(pseg_samples_)));
 
        }//end while
 
        NotifyObservers(InEGOperatingZone<EmitterType>(*this));
 
        return SuccessCode::Success;
    }

    template<typename ComplexT>
    SuccessCode ComputeCauchyApproximationOfXAtT0(Vec<ComplexT>& result)
    {
        auto& cau_times = std::get<TimeCont<ComplexT> >(cauchy_times_);
        auto& cau_samples = std::get<SampCont<ComplexT> >(cauchy_samples_);
 
        if (cau_samples.size() != this->CycleNumber() * this->EndgameSettings().num_sample_points+1)
        {
            std::stringstream err_msg;
            err_msg << "to compute cauchy approximation, cau_samples must be of size " << this->CycleNumber() * this->EndgameSettings().num_sample_points+1 << " but is of size " << cau_samples.size() << '\n';
            throw std::runtime_error(err_msg.str());
        }
 
 
        //Ensure all samples are of the same precision.
        if constexpr (tracking::TrackerTraits<TrackerType>::IsAdaptivePrec)
        {
            auto new_precision = this->EnsureAtUniformPrecision(cau_times, cau_samples);
        }
 
 
        auto total_num_pts = this->CycleNumber() * this->EndgameSettings().num_sample_points;
        auto refine_code = this->template RefineAllSamples<ComplexT>(cau_samples, cau_times);
        // Pure-(i): when the adaptive-numeric-type driver is orchestrating, a refine that double cannot
        // satisfy is a request to cross to mpfr -- propagate it so RunImplAMP migrates and retries.
        // Fixed precision (and the AMP-PowerSeries path) keep ignoring the code, exactly as before.
        if constexpr (tracking::TrackerTraits<TrackerType>::IsAdaptivePrec)
            if (this->adaptive_numeric_type_active_ &&
                (refine_code == SuccessCode::HigherPrecisionNecessary || refine_code == SuccessCode::FailedToConverge))
                return refine_code;
 
        Precision(result, Precision(cau_samples.back()));
 
        result = Vec<ComplexT>::Zero(static_cast<Eigen::Index>(this->GetSystem().NumVariables()));
        for(unsigned int ii = 0; ii < total_num_pts; ++ii)
            result += cau_samples[ii];
        result /= this->CycleNumber() * this->EndgameSettings().num_sample_points;
 
        // A NaN mean must be a FAILURE code, never Success: NaN compares false against
        // everything, so downstream convergence and security comparisons are blind to it.
        if (bertini::ContainsNaN(result))
            return SuccessCode::FailedToConverge;
 
        return SuccessCode::Success;
 
    }

    template<typename ComplexT>
    NumErrorT MinDehomNormOverLoopSamples() const
    {
        const auto& cau_samples = std::get<SampCont<ComplexT> >(cauchy_samples_);
        using std::min; using std::isnan;
        NumErrorT floor_of_loop = std::numeric_limits<NumErrorT>::max();
        for (const auto& s : cau_samples)
        {
            auto n = static_cast<NumErrorT>(this->GetSystem().InfinityNormOfDehomogenized(s));
            if (!isnan(n))
                floor_of_loop = min(floor_of_loop, n);
        }
        return floor_of_loop;
    }

    template<typename ComplexT>
    SuccessCode ComputeCauchySamples(ComplexT const& target_time)
    {
        using bertini::Precision;
 
        ClearAndSeedCauchyData<ComplexT>();
        this->CycleNumber(0);
 
 
        while( this->CycleNumber() < GetCauchySettings().fail_safe_maximum_cycle_number )
        {
            //track around the origin once.
            auto tracking_success = CircleTrack(target_time);
            this->IncrementCycleNumber(1);
 
            if(tracking_success != SuccessCode::Success)
            {
                return tracking_success;
            }
            else if(CheckClosedLoop<ComplexT>())
            {
                return SuccessCode::Success;
            }
        }
        NotifyObservers(CycleNumTooHigh<EmitterType>(*this));
        return SuccessCode::CycleNumTooHigh;
    }//end ComputeCauchySamples
 

    template<typename ComplexT>
    SuccessCode AdvanceTime(ComplexT const& target_time)
    {
        using RealT = typename Eigen::NumTraits<ComplexT>::Real;
 
        auto& ps_times = std::get<TimeCont<ComplexT> >(pseg_times_);
        auto& ps_samples = std::get<SampCont<ComplexT> >(pseg_samples_);
 
        auto& current_time = ps_times.back();
        auto& current_sample = ps_samples.back();
 
        //Generalized next_time in case if we are not trying to converge to the t = 0.
        ComplexT next_time = (target_time-current_time) * static_cast<RealT>(this->EndgameSettings().sample_factor)+current_time;
 
        if (abs(next_time - target_time) < this->EndgameSettings().min_track_time)//we are too close to t = 0 but we do not have the correct tolerance - so we exit
            return SuccessCode::MinTrackTimeReached;
 
        // advance in time
        Vec<ComplexT> next_sample;
        auto time_advance_success = this->EndgameTrackPath(next_sample,current_time, next_time, current_sample);
        if (time_advance_success != SuccessCode::Success)
        {
            NotifyObservers(EndgameFailure<EmitterType>(*this));
            return time_advance_success;
        }
 
        // Pure-(i) escalation: return BEFORE RotateOntoPS so the PSEG window stays an untouched
        // checkpoint -- RunImplAMP migrates it to mpfr and retries this advance.  Elided for fixed prec.
        if constexpr (tracking::TrackerTraits<TrackerType>::IsAdaptivePrec)
            if (this->adaptive_numeric_type_active_ &&
                this->GetTracker().GetCurrentPrecision() > this->current_endgame_precision_)
                return SuccessCode::HigherPrecisionNecessary;
 
        this->EnsureAtPrecision(next_time,Precision(next_sample));
        RotateOntoPS(next_time, next_sample);
 
        NotifyObservers(TimeAdvanced<EmitterType>(*this));
        // the sample is complete here: Cauchy does not refine the new power-series
        // sample at this site (RefineAllSamples does the window later), so the
        // event carries the sample as tracked
        // ComputedSamplePoint carries the sample/time at BaseComplexT, like CircleAdvanced.
        // The fixed/mpfr lane emits directly; the complex_dbl fast lane would have to convert
        // the point to mpfr for the event, so we only pay that when something is observing
        // (temporaries live through the synchronous NotifyObservers).
        if constexpr (std::is_same<ComplexT, BCT>::value)
            NotifyObservers(ComputedSamplePoint<EmitterType>(*this, next_sample, next_time));
        else if (this->HasObservers())
        {
            Vec<BCT> ev_pt(next_sample.size());
            for (Eigen::Index i = 0; i < next_sample.size(); ++i) ev_pt(i) = BCT(next_sample(i));
            BCT ev_t(next_time);
            NotifyObservers(ComputedSamplePoint<EmitterType>(*this, ev_pt, ev_t));
        }
        return SuccessCode::Success;
    }

    template<typename ComplexT>
    NumErrorT PoleComponentMass() const
    {
        using RealT = typename Eigen::NumTraits<ComplexT>::Real;
        auto const& samples = std::get<SampCont<ComplexT>>(cauchy_samples_);
        auto const& times   = std::get<TimeCont<ComplexT>>(cauchy_times_);
        auto const c = this->CycleNumber();
        auto const N = this->EndgameSettings().num_sample_points;
        auto const M = c * N;
        if (c == 0 || samples.size() < M + 1 || times.size() < 2)
            return NumErrorT(0);
 
        // orientation of the loop (which way theta advances); the frequency sign follows it
        RealT const orientation_test = imag(times[1] * conj(times[0]));
        RealT const orientation = (orientation_test > RealT(0)) ? RealT(1) : RealT(-1);
 
        RealT const two_pi = RealT(2) * acos(RealT(-1));
        NumErrorT mass(0);
        for (unsigned m = 1; m <= c; ++m)
        {
            Vec<ComplexT> twisted = Vec<ComplexT>::Zero(samples[0].size());
            for (unsigned j = 0; j < M; ++j)
            {
                // weight e^{+i m theta_j / c} with theta_j = orientation * 2 pi j / N:
                // picks out exactly the k = -m Puiseux mode
                RealT const angle = orientation * two_pi * RealT(j * m) / RealT(N * c);
                twisted += samples[j] * ComplexT(cos(angle), sin(angle));
            }
            twisted /= RealT(M);
            auto const twisted_norm = static_cast<NumErrorT>(twisted.template lpNorm<Eigen::Infinity>());
            if (twisted_norm > mass)
                mass = twisted_norm;
        }
        return mass;
    }

    template<typename Dummy = void>
    NumErrorT PoleComponentMassAMP() const
    {
        if (this->current_endgame_precision_ == DoublePrecision())
            return PoleComponentMass<complex_dbl>();
        return PoleComponentMass<complex_mp>();
    }

    bool PoleMassSignificant(NumErrorT mass, NumErrorT approx_scale) const
    {
        return mass > NumErrorT(1e3) * this->FinalTolerance() * (NumErrorT(1) + approx_scale);
    }

    template<typename ComplexT>
    SuccessCode RunImpl(ComplexT const& start_time, Vec<ComplexT> const& start_point, ComplexT const& target_time)
    {
        if (start_point.size()!=static_cast<Eigen::Index>(this->GetSystem().NumVariables()))
        {
            std::stringstream err_msg;
            err_msg << "number of variables in start point for CauchyEG, " << start_point.size() << ", must match the number of variables in the system, " << this->GetSystem().NumVariables();
            throw std::runtime_error(err_msg.str());
        }
 
        if (Precision(start_time)!=Precision(start_point))
        {
            std::stringstream ss;
            ss << "CauchyEG Run time and point must be of matching precision. (" << Precision(start_time) << "!=" << Precision(start_point) << ")";
            throw std::runtime_error(ss.str());
        }
 
        using RealT = typename Eigen::NumTraits<ComplexT>::Real;
 
        Vec<ComplexT>& latest_approx = this->final_approximation_;
        Vec<ComplexT>& prev_approx = this->previous_approximation_;
        NumErrorT& approx_error = this->approximate_error_;
 
 
        ClearTimesAndSamples<ComplexT>(); //clear times and samples before we begin.
        this->CycleNumber(0);
        prev_approx = start_point;
 
        auto init_success = GetIntoEGZone(start_time, start_point, target_time);
        if (init_success!= SuccessCode::Success)
            return init_success;
 
        auto cauchy_loop_success = InitialCauchyLoops<ComplexT>(target_time);
        if (cauchy_loop_success != SuccessCode::Success)
            return cauchy_loop_success;
 
 
        // The security check watches TWO honest divergence signals, taking the stronger:
        // (1) the extrapolated ENDPOINT -- a genuinely-infinite endpoint (a patched
        // projective path leaving the affine chart) has a dehomogenized norm growing to
        // infinity FASTER than any finite loop sample, the earliest signal there; and
        // (2) the loop FLOOR (min over loop samples) -- a fractional-pole diverger's
        // branches cancel around the closed loop, so its mean/endpoint stays ~0 while
        // every sample sits at ~r^(-p/q); only the floor sees it.  The floor is safe
        // where watching samples DIRECTLY was not (any single sample transiently spiking
        // above max_norm truncated finite cyclic-6 paths, 156 -> 153): a transient spike
        // lifts the max, never the min.  The pole-annihilation case -- a Laurent pole
        // whose Cauchy mean is a STATIONARY FINITE number -- is handled separately by
        // the pole-component operating-zone check below.  Tracked across consecutive
        // ARMED rounds only (see below); initialized to 0 so the check never reads
        // indeterminate values.
        RealT norm_of_dehom_prev(0), norm_of_dehom_latest(0);
 
        // Cycle-number consistency: refuse to accept a converged approximation until the cycle number
        // has reported the SAME value for num_consecutive_same_cycle_number consecutive approximations.
        // When the working precision is too low to close the loop accurately (NOT monodromy -- a
        // nonsingular endpoint is cycle 1 at every radius in high precision), the cycle number thrashes
        // (41, 14, 36, ...) and a coincidental approx_error dip must not be mistaken for convergence.
        // prev_cycle == 0 means "no prior measurement".
        unsigned prev_cycle = 0, same_cycle_count = 0;
 
        do
        {
            //Compute a cauchy approximation.  Uses the previously computed samples,
            //either from InitialCauchyLoops, or ComputeCauchySamples
            auto extrapolation_success = ComputeCauchyApproximationOfXAtT0<ComplexT>(latest_approx);
            if (extrapolation_success!=SuccessCode::Success)
                return extrapolation_success;
 
            unsigned cur_cycle = this->CycleNumber();
            if (cur_cycle == prev_cycle)
                ++same_cycle_count;
            else
                same_cycle_count = 1;
            prev_cycle = cur_cycle;
 
            approx_error = static_cast<NumErrorT>((latest_approx - prev_approx).template lpNorm<Eigen::Infinity>());
            NotifyObservers(ApproximatedRoot<EmitterType>(*this));
 
            // the operating-zone ACCEPTANCE GATE: negative-mode (pole) mass in the loop.
            // Nonzero mass means the disk between here and the target time is not clean
            // -- a pole at the target, or another branch point inside the loop -- and
            // the mean is not to be trusted, however stationary it looks.  So convergence
            // is refused until the zone is reached (pole mass insignificant).  This alone
            // cures the junk-success bug: an affine Laurent pole has a stationary FINITE
            // Cauchy mean (roots-of-unity annihilation), so no norm check can catch it,
            // but its pole mass never vanishes, so this gate never accepts it -- the path
            // instead runs to the endgame's natural terminal condition (min track time /
            // fail-safe max cycle) and returns non-Success.  There is deliberately NO
            // active pole-growth truncation: distinguishing a genuine pole (mass grows
            // unboundedly) from a branch point merely inside the loop (mass grows a round
            // or two, then dies as the radius shrinks past it) by counting growing rounds
            // produced false positives that truncated CLEAN convergent paths -- cyclic-6
            // lost genuine nonsingular solutions (156 -> 153).  The gate suffices.
            auto const pole_mass = PoleComponentMass<ComplexT>();
            bool const in_operating_zone = !PoleMassSignificant(pole_mass,
                static_cast<NumErrorT>(latest_approx.template lpNorm<Eigen::Infinity>()));
            // B1's CycleTimeCutoff semantics (manual E.5.9): below cycle_cutoff_time the
            // zone heuristics stop being consulted and the endgame behaves as if in the
            // operating zone.  Applied here to the SECURITY VALVE ONLY: a slowly-diverging
            // path (fractional-order blowup) carries significant pole mass at every radius,
            // so it is never pole-mass-in-zone and the valve never arms -- it crawls to the
            // tracker's far-larger truncation threshold (measured: a t^(-1/3) diverger took
            // ~117k extra steps over 3 more decades of |t| at up to 70 digits).  Below the
            // cutoff a genuinely convergent path's mean has long stabilized at a finite
            // value, so arming the valve cannot reintroduce the cyclic-6 false truncations.
            // Deliberately NOT applied to acceptance: forcing the zone there would re-open
            // the junk-success hole (a Laurent pole's stationary FINITE mean would be
            // accepted once below the cutoff); pole junk still runs to its honest
            // non-Success terminal.
            using std::abs;
            bool const below_cycle_cutoff =
                static_cast<NumErrorT>(abs(this->LatestTime() - target_time))
                    < GetCauchySettings().cycle_cutoff_time;
 
            // SECURITY BEFORE ACCEPTANCE, the order the power series endgame uses.  At
            // security level 0, two consecutive endpoint approximations whose dehomogenized
            // infinity norm is above max_norm truncate the path, whether or not the second
            // would also have been accepted.  Asked the other way round, a path that
            // converged to a finite point above the ceiling came back Success here and
            // SecurityMaxNormReached from the power series endgame -- the two endgames must
            // differ in how they estimate the root and in nothing else.
            //
            // Security truncation happens ONLY in the operating zone.  Outside it the
            // loop may encircle a pole or branch point and the Cauchy mean is garbage --
            // its norm means nothing, and truncating on it killed clean convergent paths
            // (the cyclic-6 residual, 156 -> 155).  Inside the zone the mean is
            // trustworthy: a genuinely-infinite endpoint is in-zone AND large (clean
            // homogeneous convergence, the endpoint growing to infinity faster than any
            // sample), so two consecutive IN-ZONE rounds above max_norm truncate honestly.
            // An out-of-zone round restarts the two-consecutive count.
            if (this->SecuritySettings().level <= 0)
            {
                if (in_operating_zone || below_cycle_cutoff)
                {
                    // Watch the STRONGER of two honest divergence signals: the endpoint (the
                    // earliest signal for a clean infinite endpoint) and the loop FLOOR (min
                    // over loop samples) -- the mean is structurally blind to fractional-pole
                    // divergers, whose branches cancel around the closed loop while every
                    // sample sits at ~r^(-p/q).  See MinDehomNormOverLoopSamples.
                    using std::max;
                    norm_of_dehom_latest = max(
                        this->GetSystem().InfinityNormOfDehomogenized(latest_approx),
                        static_cast<RealT>(MinDehomNormOverLoopSamples<ComplexT>()));
                    if (this->BeyondSecurityMaxNorm(norm_of_dehom_prev) &&
                        this->BeyondSecurityMaxNorm(norm_of_dehom_latest))
                    {
                        NotifyObservers(SecurityMaxNormReached<EmitterType>(*this));
                        return SuccessCode::SecurityMaxNormReached;
                    }
                    norm_of_dehom_prev = norm_of_dehom_latest;
                }
                else
                    norm_of_dehom_prev = RealT(0);   // out of zone: restart the count
            }
 
            // An approximation above the ceiling is never ACCEPTED at security level 0: it is
            // the first of the two that truncate, or the second.  Without this the count above
            // (armed rounds only) could still be at one when acceptance, which needs a single
            // in-zone round, let the path through as a Success above max_norm.
            bool const above_ceiling = this->SecuritySettings().level <= 0
                && this->BeyondSecurityMaxNorm(
                        this->GetSystem().InfinityNormOfDehomogenized(latest_approx));
 
            if (in_operating_zone
                && !above_ceiling
                && approx_error < this->FinalTolerance()
                && same_cycle_count >= GetCauchySettings().num_consecutive_same_cycle_number)
            {
                NotifyObservers(Converged<EmitterType>(*this));
                return SuccessCode::Success;
            }
 
            prev_approx = latest_approx;
 
            auto advance_success = AdvanceTime<ComplexT>(target_time);
            if (advance_success != SuccessCode::Success)
                return advance_success;
 
            // then compute the next set of cauchy samples used for extrapolating the point at target time
            auto cauchy_samples_success = ComputeCauchySamples(target_time);
            if (cauchy_samples_success != SuccessCode::Success)
                return cauchy_samples_success;
 
        } while (true);
 
        return SuccessCode::Success;
    } //end main CauchyEG function
 
 
    // ================================================================================================
    //   Adaptive-numeric-type (double-first) Cauchy endgame.
    //
    //   Computes in the hardware-complex_dbl fast lane while the AMP tracker's authoritative precision
    //   stays double, and crosses to complex_mp only when the tracker escalates a TrackPath/Refine
    //   (pure-(i)).  Fixed precision never enters here -- base Run() sends it to RunImpl<BCT>.  Every
    //   method below is a member template so the explicit fixed-precision class instantiations do not
    //   force-compile them (which would std::get a complex_mp/complex_dbl slot a single-type endgame
    //   does not have).  The complex_dbl->complex_mp container migration is the one piece with no
    //   pre-existing analog; mpfr->higher-mpfr co-vary already happens via EnsureAtUniformPrecision.
    // ================================================================================================

    template<typename Dummy = void>
    SuccessCode RunImplAMP(BCT const& start_time, Vec<BCT> const& start_point, BCT const& target_time)
    {
        using bertini::Precision;
        using RealT = typename Eigen::NumTraits<BCT>::Real;
 
        if (start_point.size()!=static_cast<Eigen::Index>(this->GetSystem().NumVariables()))
        {
            std::stringstream err_msg;
            err_msg << "number of variables in start point for CauchyEG, " << start_point.size() << ", must match the number of variables in the system, " << this->GetSystem().NumVariables();
            throw std::runtime_error(err_msg.str());
        }
 
        // Arm the escalation hooks in the shared phase methods (CircleTrack / AdvanceTime / refine);
        // disarm on every exit path.
        this->adaptive_numeric_type_active_ = true;
        struct Disarmer { bool& flag; ~Disarmer(){ flag = false; } } disarm{this->adaptive_numeric_type_active_};
 
        // Both lanes start empty: a phase clears only the lane it runs in, so a path whose endgame
        // stays in double would otherwise leave an earlier path's multiprecision samples and times
        // in place for anything that reads that lane first (#378).
        ClearTimesAndSamples<complex_dbl>();
        ClearTimesAndSamples<complex_mp>();
 
        // Start in the precision the tracker handed us at the endgame boundary: double for the easy
        // majority, already-mpfr for the few paths that escalated before the endgame.  Migrate UP only.
        this->current_endgame_precision_ = std::max(DoublePrecision(), Precision(start_point));
        this->CycleNumber(0);
 
        // ---- INITIALIZATION (GetIntoEGZone + InitialCauchyLoops).  On escalation, restart init from the
        //      boundary at the higher precision: init lives at the well-conditioned large-|t| end where
        //      escalation is rare, and a fresh init clears its own containers so it is self-consistent. ----
        while (true)
        {
            SuccessCode init_code;
            if (this->current_endgame_precision_ == DoublePrecision())
                init_code = RunInitSegmentT<complex_dbl>(complex_dbl(start_time), this->DowncastToDouble(start_point), complex_dbl(target_time));
            else
                init_code = RunInitSegmentT<complex_mp>(this->AtActivePrecisionScalar(start_time), this->AtActivePrecisionVec(start_point), this->AtActivePrecisionScalar(target_time));
 
            if (init_code == SuccessCode::HigherPrecisionNecessary)
            {
                auto const previous_precision = this->current_endgame_precision_;
                this->current_endgame_precision_ = this->NextEscalatedPrecision();
                SetThreadPrecision(this->current_endgame_precision_);
                // the samples announced at the lower precision are superseded (the window is
                // tracked again at the new one); say so, so a collector drops them as we do
                NotifyObservers(PrecisionChanged<EmitterType>(*this, previous_precision, this->current_endgame_precision_));
                NotifyObservers(SamplesRecomputedAtHigherPrecision<EmitterType>(*this));
                continue;
            }
            if (init_code != SuccessCode::Success)
                return init_code;
            break;
        }
 
        // ---- MAIN CONVERGENCE LOOP.  Each phase self-heals: on escalation it migrates the durable state
        //      up to mpfr (widening the retained PSEG window, never re-tracking it) and retries in mpfr.
        //      The approximations live in BCT, so the per-iteration bookkeeping arithmetic is mpfr -- one
        //      small vector op next to the tracking, which itself stays in the fast lane. ----
        // the security check watches the extrapolated ENDPOINT for honest divergence to
        // infinity (which grows faster than any finite loop sample), but ONLY in the
        // operating zone (a trustworthy mean) -- see RunImpl for the full reasoning.
        // Acceptance is gated on the same operating-zone measurement (the finite-mean
        // pole case).  Tracked across consecutive in-zone rounds; starts at 0.
        RealT norm_of_dehom_prev(0), norm_of_dehom_latest(0);
 
        unsigned prev_cycle = 0, same_cycle_count = 0;
 
        while (true)
        {
            auto extrap_code = ComputeCauchyApproxAMP();
            if (extrap_code != SuccessCode::Success)
                return extrap_code;
 
            unsigned cur_cycle = this->CycleNumber();
            if (cur_cycle == prev_cycle) ++same_cycle_count; else same_cycle_count = 1;
            prev_cycle = cur_cycle;
 
            Precision(this->previous_approximation_, Precision(this->final_approximation_));
            this->approximate_error_ = static_cast<NumErrorT>((this->final_approximation_ - this->previous_approximation_).template lpNorm<Eigen::Infinity>());
            NotifyObservers(ApproximatedRoot<EmitterType>(*this));
 
            // the operating-zone ACCEPTANCE GATE (see RunImpl for the full reasoning):
            // convergence is refused while pole mass is significant, which alone cures
            // junk-success; there is deliberately no active pole-growth truncation (it
            // false-positived on clean convergent paths -- cyclic-6, 156 -> 153).
            auto const pole_mass = PoleComponentMassAMP();
            bool const in_operating_zone = !PoleMassSignificant(pole_mass,
                static_cast<NumErrorT>(this->final_approximation_.template lpNorm<Eigen::Infinity>()));
            // B1's CycleTimeCutoff semantics (manual E.5.9): below cycle_cutoff_time the
            // zone heuristics stop being consulted and the endgame behaves as if in the
            // operating zone.  Applied here to the SECURITY VALVE ONLY: a slowly-diverging
            // path (fractional-order blowup) carries significant pole mass at every radius,
            // so it is never pole-mass-in-zone and the valve never arms -- it crawls to the
            // tracker's far-larger truncation threshold (measured: a t^(-1/3) diverger took
            // ~117k extra steps over 3 more decades of |t| at up to 70 digits).  Below the
            // cutoff a genuinely convergent path's mean has long stabilized at a finite
            // value, so arming the valve cannot reintroduce the cyclic-6 false truncations.
            // Deliberately NOT applied to acceptance: forcing the zone there would re-open
            // the junk-success hole (a Laurent pole's stationary FINITE mean would be
            // accepted once below the cutoff); pole junk still runs to its honest
            // non-Success terminal.
            using std::abs;
            bool const below_cycle_cutoff =
                static_cast<NumErrorT>(abs(this->LatestTime() - this->AtActivePrecisionScalar(target_time)))
                    < GetCauchySettings().cycle_cutoff_time;
 
            // security BEFORE acceptance, as in the power series endgame; truncation ONLY in
            // the operating zone (a trustworthy mean); an out-of-zone round restarts the
            // two-consecutive count.  See RunImpl.
            if (this->SecuritySettings().level <= 0)
            {
                if (in_operating_zone || below_cycle_cutoff)
                {
                    // The stronger of endpoint norm and loop floor; the loop samples live in
                    // whichever lane is active.  See RunImpl and MinDehomNormOverLoopSamples.
                    NumErrorT loop_floor = (this->current_endgame_precision_ == DoublePrecision())
                        ? MinDehomNormOverLoopSamples<complex_dbl>()
                        : MinDehomNormOverLoopSamples<complex_mp>();
                    using std::max;
                    norm_of_dehom_latest = max(
                        this->GetSystem().InfinityNormOfDehomogenized(this->final_approximation_),
                        static_cast<RealT>(loop_floor));
                    if (this->BeyondSecurityMaxNorm(norm_of_dehom_prev) &&
                        this->BeyondSecurityMaxNorm(norm_of_dehom_latest))
                    {
                        NotifyObservers(SecurityMaxNormReached<EmitterType>(*this));
                        return SuccessCode::SecurityMaxNormReached;
                    }
                    norm_of_dehom_prev = norm_of_dehom_latest;
                }
                else
                    norm_of_dehom_prev = RealT(0);   // out of zone: restart the count
            }
 
            // never accept an approximation above the ceiling at security level 0 (see RunImpl)
            bool const above_ceiling = this->SecuritySettings().level <= 0
                && this->BeyondSecurityMaxNorm(
                        this->GetSystem().InfinityNormOfDehomogenized(this->final_approximation_));
 
            if (in_operating_zone
                && !above_ceiling
                && this->approximate_error_ < this->FinalTolerance()
                && same_cycle_count >= GetCauchySettings().num_consecutive_same_cycle_number)
            {
                NotifyObservers(Converged<EmitterType>(*this));
                return SuccessCode::Success;
            }
 
            this->previous_approximation_ = this->final_approximation_;
 
            auto advance_code = AdvanceTimeAMP(target_time);
            if (advance_code != SuccessCode::Success)
                return advance_code;
 
            auto samples_code = ComputeCauchySamplesAMP(target_time);
            if (samples_code != SuccessCode::Success)
                return samples_code;
        }
 
        return SuccessCode::Success;
    }
 
 
    // Initialization at one numeric type: seed the PSEG window and reach the EG operating zone, exactly
    // as the head of RunImpl does, but reporting HigherPrecisionNecessary up to the driver on escalation.
    template<typename ComplexT>
    SuccessCode RunInitSegmentT(ComplexT const& start_time, Vec<ComplexT> const& start_point, ComplexT const& target_time)
    {
        ClearTimesAndSamples<ComplexT>();
        this->CycleNumber(0);
        this->previous_approximation_ = this->ToBCT(start_point, this->current_endgame_precision_);
 
        auto init_success = GetIntoEGZone(start_time, start_point, target_time);
        if (init_success != SuccessCode::Success)
            return init_success;
 
        return InitialCauchyLoops<ComplexT>(target_time);
    }
 
 
    // Cauchy mean (extrapolation) at the active numeric type, written into final_approximation_ (BCT).
    template<typename Dummy = void>
    SuccessCode ComputeCauchyApproxAMP()
    {
        unsigned guard = 0;
        while (true)
        {
            SuccessCode code;
            if (this->current_endgame_precision_ == DoublePrecision())
            {
                Vec<complex_dbl> r;
                code = ComputeCauchyApproximationOfXAtT0<complex_dbl>(r);
                if (code == SuccessCode::Success)
                    this->final_approximation_ = this->ToBCT(r, this->current_endgame_precision_);
            }
            else
            {
                code = ComputeCauchyApproximationOfXAtT0<complex_mp>(this->final_approximation_);
            }
 
            if (code == SuccessCode::HigherPrecisionNecessary || code == SuccessCode::FailedToConverge)
            {
                if (this->template EscalateAndMigrate<>(++guard) != SuccessCode::Success)
                    return SuccessCode::HigherPrecisionNecessary;
                continue;
            }
            return code;
        }
    }
 
 
    // Advance time at the active numeric type; migrate-and-retry in mpfr on escalation.  AdvanceTime
    // returns HigherPrecisionNecessary BEFORE it rotates the PSEG window, so the window stays a clean
    // checkpoint and the retry continues from the widened window with no double-advance.
    template<typename Dummy = void>
    SuccessCode AdvanceTimeAMP(BCT const& target_time)
    {
        unsigned guard = 0;
        while (true)
        {
            SuccessCode code = (this->current_endgame_precision_ == DoublePrecision())
                ? AdvanceTime<complex_dbl>(complex_dbl(target_time))
                : AdvanceTime<complex_mp>(this->AtActivePrecisionScalar(target_time));
 
            if (code == SuccessCode::HigherPrecisionNecessary)
            {
                if (this->template EscalateAndMigrate<>(++guard) != SuccessCode::Success)
                    return SuccessCode::HigherPrecisionNecessary;
                continue;
            }
            return code;
        }
    }
 
 
    // Build a closed Cauchy loop's samples at the active numeric type; migrate-and-retry in mpfr on
    // escalation.  ComputeCauchySamples clears and re-seeds its cauchy data from the PSEG window, so the
    // mpfr retry simply re-tracks the circle from the (migrated) window -- it never reuses the lossy,
    // double-tracked partial circle that triggered the escalation.
    template<typename Dummy = void>
    SuccessCode ComputeCauchySamplesAMP(BCT const& target_time)
    {
        unsigned guard = 0;
        while (true)
        {
            SuccessCode code = (this->current_endgame_precision_ == DoublePrecision())
                ? ComputeCauchySamples<complex_dbl>(complex_dbl(target_time))
                : ComputeCauchySamples<complex_mp>(this->AtActivePrecisionScalar(target_time));
 
            if (code == SuccessCode::HigherPrecisionNecessary)
            {
                if (this->template EscalateAndMigrate<>(++guard) != SuccessCode::Success)
                    return SuccessCode::HigherPrecisionNecessary;
                continue;
            }
            return code;
        }
    }
 
 
    // Flavor-specific: cross every Cauchy container from the complex_dbl slot to complex_mp (or, if
    // already mpfr, raise its precision uniformly), via the shared base Cross* / SetPrecision helpers.
    // Widen-only by default -- retained samples were already tracked/refined to final_tolerance
    // (pure-(i)/B).  Called by the base EscalateAndMigrate.
    template<typename Dummy = void>
    void MigrateContainersToPrecision(unsigned newprec)
    {
        using bertini::Precision;
        if (this->current_endgame_precision_ == DoublePrecision())
        {
            this->template CrossTimesUp<>(pseg_times_,     newprec);
            this->template CrossSampsUp<>(pseg_samples_,   newprec);
            this->template CrossTimesUp<>(cauchy_times_,   newprec);
            this->template CrossSampsUp<>(cauchy_samples_, newprec);
        }
        else
        {
            tracking::adaptive::SetPrecision(std::get<TimeCont<complex_mp>>(pseg_times_),     newprec);
            tracking::adaptive::SetPrecision(std::get<SampCont<complex_mp>>(pseg_samples_),   newprec);
            tracking::adaptive::SetPrecision(std::get<TimeCont<complex_mp>>(cauchy_times_),   newprec);
            tracking::adaptive::SetPrecision(std::get<SampCont<complex_mp>>(cauchy_samples_), newprec);
        }
 
        if (this->final_approximation_.size()    > 0) Precision(this->final_approximation_,    newprec);
        if (this->previous_approximation_.size() > 0) Precision(this->previous_approximation_, newprec);
 
        if (this->EndgameSettings().refine_when_increasing_precision)
        {
            auto& cau_t = std::get<TimeCont<complex_mp>>(cauchy_times_);
            auto& cau_s = std::get<SampCont<complex_mp>>(cauchy_samples_);
            if (!cau_s.empty())
                this->template RefineAllSamples<complex_mp>(cau_s, cau_t);
        }
    }
 
};
 
 
}} // namespaces