config.hpp

Go to the documentation of this file.

//This file is part of Bertini 2.
//
//trackers/include/bertini2/trackers/config.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.
//
//trackers/include/bertini2/trackers/config.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 tracking/include/bertini2/trackers/config.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.
 
#ifndef BERTINI_TRACKING_CONFIG_HPP
#define BERTINI_TRACKING_CONFIG_HPP

#include "bertini2/mpfr_extensions.hpp"
#include "bertini2/eigen_extensions.hpp"
#include "bertini2/system/system.hpp"
#include "bertini2/detail/typelist.hpp"
 
#include "bertini2/common/config.hpp"
 
 
namespace bertini
{
 
namespace tracking{
 

    enum class PrecisionType //E.2.1
    {
        Fixed,
        FixedMultiple,
        Adaptive
    };
 

    enum class Predictor //E.4.3
    {
        Constant,
        Euler,
        Heun,
        RK4,
        HeunEuler,
        RKNorsett34,
        RKF45,
        RKCashKarp45,
        RKDormandPrince56,
        RKVerner67
    };
 
 
 
 
 
 

    struct StepMetadata
    {
        NumErrorT norm_J = 0;                    
        NumErrorT norm_J_inverse = 0;            
        NumErrorT condition_number_estimate = 0; 
        NumErrorT size_proportion = 0;           
        NumErrorT error_estimate = 0;            
        NumErrorT norm_delta_z = 0;              
    };
 

    struct SteppingConfig
    {
        // mpq_rational: exact rationals with no MPFR precision state — safe in DefaultConstruct<T>::value statics.
        // real_mp fields here would be initialized at BMP's startup precision (20) and contaminate
        // tracker arithmetic when target precision < 20 via preserve_related_precision.
        mpq_rational initial_step_size{1, 10}; 
        mpq_rational max_step_size{1, 10};     
        double       min_step_size = 1e-100;   
 
        mpq_rational step_size_success_factor{2, 1}; 
        mpq_rational step_size_fail_factor{1, 2};    
 
        unsigned consecutive_successful_steps_before_stepsize_increase = 5; 
 
        unsigned min_num_steps = 1; 
        unsigned max_num_steps = 1e5; 
 
        unsigned frequency_of_CN_estimation = 1; 
    };
 
 

    struct NewtonConfig
    {
        unsigned max_num_newton_iterations = 2; 
        unsigned min_num_newton_iterations = 1; 
    };
 

    struct TrackerConfig
    {
        Predictor predictor = Predictor::RKF45; 
        NumErrorT tracking_tolerance = NumErrorT(1e-5); 
        NumErrorT path_truncation_threshold = NumErrorT(1e5); 
    };
 
 
 
 
 
 
 

    struct FixedPrecisionConfig
    {
        using RealT = double;  

        unsigned precision = 0;

        explicit
        FixedPrecisionConfig(System const& /*sys*/)
        { }
 
        FixedPrecisionConfig() = default;
    };
 

    inline
    std::ostream& operator<<(std::ostream & out, FixedPrecisionConfig const& /*fpc*/)
    {
        return out;
    }
 

    struct AdaptiveMultiplePrecisionConfig
    {
        NumErrorT coefficient_bound;  
        NumErrorT degree_bound; 

        NumErrorT linear_solve_error_bound;

        NumErrorT jacobian_eval_error_bound;

        NumErrorT function_eval_error_bound;
 
        int safety_digits_1 = 1; 
        int safety_digits_2 = 1; 
        unsigned int maximum_precision = 300; 
 
        // Note: a single setting -- Bertini 1's StepsForIncrease, i.e.
        // SteppingConfig::consecutive_successful_steps_before_stepsize_increase -- gates BOTH stepsize
        // increase AND precision decrease (the required number of consecutive successful steps).
        // Precision decrease is additionally subject to B1's digits-margin hysteresis (see
        // ExtraDigitsBeforePrecisionDecrease in amp_tracker.hpp).  The old, duplicate AMP-config setting
        // `consecutive_successful_steps_before_precision_decrease` (a B2 deviation) was removed.
 
        unsigned max_num_precision_decreases = 10; 
 

        void SetBoundsFrom(System const& sys)
        {
            using std::pow;
 
            // Refuse before DegreeBound() does, so the message is about what the caller was
            // trying to do rather than about degrees.  The two bounds are on the error of
            // evaluating the Jacobian and the functions; the degree and coefficient bounds are
            // merely the recipe for them that holds for polynomials.  There is no accepted
            // recipe for an analytic system -- see issue #439 -- so we do not invent one.
            if (!sys.IsPolynomial())
                throw std::runtime_error("adaptive precision derives its error bounds from the "
                    "degree of the system, and this system is not a polynomial one, so it has no "
                    "degree.  Two ways on: track with a fixed-precision tracker, at double or at "
                    "multiple precision, which needs no such bound; or set this config's "
                    "jacobian_eval_error_bound and function_eval_error_bound yourself and hand it "
                    "to the tracker, choosing values you can defend for your system.");
 
            linear_solve_error_bound = pow(NumErrorT(sys.NumVariables()),2);
            degree_bound = sys.DegreeBound();
            coefficient_bound = sys.CoefficientBound<complex_dbl>();
        }
 

        void SetErrorBoundsFromDegreeAndCoefficient()
        {
            jacobian_eval_error_bound = degree_bound*(degree_bound-NumErrorT(1))*coefficient_bound;
            function_eval_error_bound = degree_bound*coefficient_bound;
        }

        void SetAMPConfigFrom(System const& sys)
        {
            SetBoundsFrom(sys);
            SetErrorBoundsFromDegreeAndCoefficient();
        }

        AdaptiveMultiplePrecisionConfig() : coefficient_bound(1000), degree_bound(5), linear_solve_error_bound(1), safety_digits_1(1), safety_digits_2(1), maximum_precision(300)
        {
            SetErrorBoundsFromDegreeAndCoefficient();
        }

        explicit
        AdaptiveMultiplePrecisionConfig(System const& sys) : AdaptiveMultiplePrecisionConfig()
        {
            if (sys.IsPolynomial())
                SetAMPConfigFrom(sys);
        }
    }; // re: AdaptiveMultiplePrecisionConfig

    inline
    std::ostream& operator<<(std::ostream & out, AdaptiveMultiplePrecisionConfig const& AMP)
    {
        out << "coefficient_bound: " << AMP.coefficient_bound << "\n";
        out << "degree_bound: " << AMP.degree_bound << "\n";
        out << "linear_solve_error_bound: " << AMP.linear_solve_error_bound << "\n";
        out << "jacobian_eval_error_bound: " << AMP.jacobian_eval_error_bound << "\n";
        out << "function_eval_error_bound: " << AMP.function_eval_error_bound << "\n";
        out << "safety_digits_1: " << AMP.safety_digits_1 << "\n";
        out << "safety_digits_2: " << AMP.safety_digits_2 << "\n";
        out << "max_num_precision_decreases: " << AMP.max_num_precision_decreases << "\n";
        return out;
    }
 

    inline
    static
    AdaptiveMultiplePrecisionConfig AMPConfigFrom(System const& sys)
    {
        AdaptiveMultiplePrecisionConfig AMP;
        AMP.SetAMPConfigFrom(sys);
        return AMP;
    }
 
    // forward declarations
    template<class D>
    class Tracker;
    template<class D>
    class FixedPrecisionTracker;
    class MultiplePrecisionTracker;
    class DoublePrecisionTracker;
 
 
// now for the TrackerTraits structs, which enable lookup of correct settings objects and types, etc.
    template<class T>
    struct TrackerTraits
    {};
 

 
    template<>
    struct TrackerTraits<DoublePrecisionTracker>
    {
        static constexpr char const* kRecordName = "double";  
        using BaseComplexT = complex_dbl;
        using BaseRealT = double;
        using EventEmitterType = FixedPrecisionTracker<DoublePrecisionTracker>;
        using PrecisionConfig = FixedPrecisionConfig;
        enum {
            IsFixedPrec = 1,
            IsAdaptivePrec = 0
        };
 
        using NeededTypes = detail::TypeList<complex_dbl>;
        using NeededConfigs = detail::TypeList<
            SteppingConfig,
            NewtonConfig,
            TrackerConfig,
            PrecisionConfig
            >;
    };
 
 
    template<>
    struct TrackerTraits<MultiplePrecisionTracker>
    {
        static constexpr char const* kRecordName = "multiple";  
        using BaseComplexT = complex_mp;
        using BaseRealT = real_mp;
        using EventEmitterType = FixedPrecisionTracker<MultiplePrecisionTracker>;
        using PrecisionConfig = FixedPrecisionConfig;
 
        enum {
            IsFixedPrec = 1,
            IsAdaptivePrec = 0
        };
 
        using NeededTypes = detail::TypeList<complex_mp>;
 
        using NeededConfigs = detail::TypeList<
            SteppingConfig,
            NewtonConfig,
            TrackerConfig,
            PrecisionConfig
            >;
    };
 
 
    class AMPTracker; // forward declare
    template<>
    struct TrackerTraits<AMPTracker>
    {
        static constexpr char const* kRecordName = "adaptive";  
        using BaseComplexT = complex_mp;
        using BaseRealT = real_mp;
        using EventEmitterType = AMPTracker;
        using PrecisionConfig = AdaptiveMultiplePrecisionConfig;
 
        enum {
            IsFixedPrec = 0,
            IsAdaptivePrec = 1
        };
 
        using NeededTypes = detail::TypeList<complex_dbl, complex_mp>;
 
        using NeededConfigs = detail::TypeList<
            SteppingConfig,
            NewtonConfig,
            TrackerConfig,
            PrecisionConfig
            >;
    };
 
 
 
 
    template<class D>
    struct TrackerTraits<FixedPrecisionTracker<D> > : public TrackerTraits<D>
    {
        using BaseComplexT = typename TrackerTraits<D>::BaseComplexT;
        using BaseRealT = typename TrackerTraits<D>::BaseRealT;
        using EventEmitterType = typename TrackerTraits<D>::EventEmitterType;
        using PrecisionConfig = typename TrackerTraits<D>::PrecisionConfig;
 
        enum {
            IsFixedPrec = 0,
            IsAdaptivePrec = 1
        };
 
        using NeededTypes = typename TrackerTraits<D>::NeededTypes;
        using NeededConfigs = typename TrackerTraits<D>::NeededConfigs;
    };

 
} // re: namespace tracking
} // re: namespace bertini
 
 
#endif