854 lines
24 KiB
C++
854 lines
24 KiB
C++
//This file is part of Bertini 2.
|
|
//
|
|
//euler_test.cpp 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.
|
|
//
|
|
//euler_test.cpp 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 euler_test.cpp. If not, see <http://www.gnu.org/licenses/>.
|
|
//
|
|
// Copyright(C) 2015 - 2021 by 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.
|
|
|
|
// individual authors of this file include:
|
|
// silviana amethyst, university of wisconsin eau claire
|
|
|
|
|
|
|
|
|
|
|
|
#include <boost/test/unit_test.hpp>
|
|
|
|
#include <boost/multiprecision/mpfr.hpp>
|
|
#include "bertini2/mpfr_complex.hpp"
|
|
|
|
#include "bertini2/trackers/ode_predictors.hpp"
|
|
|
|
|
|
|
|
|
|
extern double threshold_clearance_d;
|
|
extern bertini::mpfr_float threshold_clearance_mp;
|
|
extern unsigned TRACKING_TEST_MPFR_DEFAULT_DIGITS;
|
|
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_SUITE(euler_predict_tracking_basics)
|
|
|
|
using System = bertini::System;
|
|
using Variable = bertini::node::Variable;
|
|
using Float = bertini::node::Float;
|
|
using ExplicitRKPredictor = bertini::tracking::predict::ExplicitRKPredictor;
|
|
|
|
using Var = std::shared_ptr<Variable>;
|
|
|
|
using VariableGroup = bertini::VariableGroup;
|
|
|
|
|
|
using dbl = std::complex<double>;
|
|
using mpfr = bertini::mpfr_complex;
|
|
using mpfr_float = bertini::mpfr_float;
|
|
|
|
|
|
template<typename NumType> using Vec = bertini::Vec<NumType>;
|
|
template<typename NumType> using Mat = bertini::Mat<NumType>;
|
|
using bertini::DefaultPrecision;
|
|
using bertini::Precision;
|
|
|
|
BOOST_AUTO_TEST_CASE(circle_line_euler_double)
|
|
{
|
|
|
|
|
|
// Starting point in spacetime step
|
|
Vec<dbl> current_space(2);
|
|
current_space << dbl(2.3,0.2), dbl(1.1, 1.87);
|
|
|
|
// Starting time
|
|
dbl current_time(0.9);
|
|
// Time step
|
|
dbl delta_t(-0.1);
|
|
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x + 5*y) );
|
|
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
BOOST_CHECK_EQUAL(AMP.degree_bound,2);
|
|
AMP.coefficient_bound = 5;
|
|
|
|
double norm_J, norm_J_inverse, size_proportion;
|
|
|
|
Vec<dbl> predicted(2);
|
|
predicted << dbl(2.40310963516214640018253210912048,0.187706567388887830930493342816564),
|
|
dbl(0.370984337833979085688209698697074, 1.30889906180158745272421523674049);
|
|
|
|
Vec<dbl> euler_prediction_result;
|
|
|
|
double tracking_tolerance(1e-5);
|
|
double condition_number_estimate;
|
|
unsigned num_steps_since_last_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(euler_prediction_result,
|
|
size_proportion,
|
|
norm_J, norm_J_inverse,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_condition_number_computation,
|
|
frequency_of_CN_estimation,
|
|
tracking_tolerance,
|
|
AMP);
|
|
|
|
BOOST_CHECK(success_code==bertini::SuccessCode::Success);
|
|
BOOST_CHECK_EQUAL(euler_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < euler_prediction_result.size(); ++ii)
|
|
BOOST_CHECK(abs(euler_prediction_result(ii)-predicted(ii)) < threshold_clearance_d);
|
|
|
|
}
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(circle_line_euler_mp)
|
|
{
|
|
bertini::DefaultPrecision(TRACKING_TEST_MPFR_DEFAULT_DIGITS);
|
|
|
|
// Starting point in spacetime step
|
|
Vec<mpfr> current_space(2);
|
|
current_space << mpfr("2.3","0.2"), mpfr("1.1", "1.87");
|
|
|
|
// Starting time
|
|
mpfr current_time("0.9");
|
|
// Time step
|
|
mpfr delta_t("-0.1");
|
|
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x + 5*y) );
|
|
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
BOOST_CHECK_EQUAL(AMP.degree_bound,2);
|
|
AMP.coefficient_bound = 5;
|
|
|
|
double norm_J, norm_J_inverse, size_proportion;
|
|
|
|
Vec<mpfr> predicted(2);
|
|
predicted << mpfr("2.40310963516214640018253210912048","0.187706567388887830930493342816564"),
|
|
mpfr("0.370984337833979085688209698697074", "1.30889906180158745272421523674049");
|
|
|
|
Vec<mpfr> euler_prediction_result;
|
|
|
|
double tracking_tolerance = 1e-5;
|
|
double condition_number_estimate;
|
|
unsigned num_steps_since_last_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(euler_prediction_result,
|
|
size_proportion,
|
|
norm_J, norm_J_inverse,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_condition_number_computation,
|
|
frequency_of_CN_estimation,
|
|
tracking_tolerance,
|
|
AMP);
|
|
|
|
BOOST_CHECK(success_code==bertini::SuccessCode::Success);
|
|
BOOST_CHECK_EQUAL(euler_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < euler_prediction_result.size(); ++ii)
|
|
BOOST_CHECK(abs(euler_prediction_result(ii)-predicted(ii)) < threshold_clearance_mp);
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(monodromy_euler_d)
|
|
{
|
|
bertini::DefaultPrecision(TRACKING_TEST_MPFR_DEFAULT_DIGITS);
|
|
|
|
// Starting point in spacetime step
|
|
Vec<dbl> current_space(2);
|
|
current_space << dbl(4.641588833612776e-1), dbl(7.416198487095662e-1);
|
|
|
|
// Starting time
|
|
dbl current_time(0.7);
|
|
// Time step
|
|
dbl delta_t(-0.01);
|
|
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,3)-1) + (1-t)*(pow(x,3) + 2) );
|
|
sys.AddFunction( t*(pow(y,2)-1) + (1-t)*(pow(y,2) + mpfr_float("0.5")) );
|
|
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
BOOST_CHECK_EQUAL(AMP.degree_bound,3);
|
|
AMP.coefficient_bound = 2;
|
|
|
|
|
|
Vec<dbl> predicted(2);
|
|
predicted << dbl(0.417742995025149735840732480384),
|
|
dbl(0.731506850772617663577442383933525);
|
|
|
|
Vec<dbl> euler_prediction_result;
|
|
|
|
double tracking_tolerance(1e-5);
|
|
double condition_number_estimate;
|
|
unsigned num_steps_since_last_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(euler_prediction_result,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_condition_number_computation,
|
|
frequency_of_CN_estimation,
|
|
tracking_tolerance);
|
|
|
|
BOOST_CHECK(success_code==bertini::SuccessCode::Success);
|
|
BOOST_CHECK_EQUAL(euler_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < euler_prediction_result.size(); ++ii)
|
|
BOOST_CHECK(abs(euler_prediction_result(ii)-predicted(ii)) < threshold_clearance_d);
|
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(monodromy_euler_mp)
|
|
{
|
|
bertini::DefaultPrecision(TRACKING_TEST_MPFR_DEFAULT_DIGITS);
|
|
|
|
// Starting point in spacetime step
|
|
Vec<mpfr> current_space(2);
|
|
current_space << mpfr("0.464158883361277585510862309093"), mpfr("0.74161984870956629487113974408");
|
|
|
|
// Starting time
|
|
mpfr current_time("0.7");
|
|
// Time step
|
|
mpfr delta_t("-0.01");
|
|
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,3)-1) + (1-t)*(pow(x,3) + 2) );
|
|
sys.AddFunction( t*(pow(y,2)-1) + (1-t)*(pow(y,2) + mpfr_float("0.5")) );
|
|
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
BOOST_CHECK_EQUAL(AMP.degree_bound,3);
|
|
AMP.coefficient_bound = 2;
|
|
|
|
|
|
Vec<mpfr> predicted(2);
|
|
predicted << mpfr("0.417742995025149735840732480384"),
|
|
mpfr("0.731506850772617663577442383934");
|
|
|
|
Vec<mpfr> euler_prediction_result;
|
|
|
|
double tracking_tolerance = 1e-5;
|
|
double condition_number_estimate;
|
|
unsigned num_steps_since_last_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(euler_prediction_result,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_condition_number_computation,
|
|
frequency_of_CN_estimation,
|
|
tracking_tolerance);
|
|
|
|
BOOST_CHECK(success_code==bertini::SuccessCode::Success);
|
|
BOOST_CHECK_EQUAL(euler_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < euler_prediction_result.size(); ++ii)
|
|
BOOST_CHECK(abs(euler_prediction_result(ii)-predicted(ii)) < threshold_clearance_mp);
|
|
}
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(euler_predict_linear_algebra_fails_d)
|
|
{
|
|
// Circle line homotopy has singular point at (x,y) = (1,-4) and t = .75
|
|
|
|
// Starting point in spacetime step
|
|
Vec<dbl> current_space(2);
|
|
current_space << dbl(1.0), dbl(-4.0);
|
|
|
|
// Starting time
|
|
dbl current_time(.75+1e-13);
|
|
// Time step
|
|
dbl delta_t(-0.1);
|
|
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x - 5*y) );
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
AMP.coefficient_bound = 5;
|
|
|
|
double tracking_tolerance(1e-5);
|
|
double condition_number_estimate;
|
|
|
|
unsigned num_steps_since_last_cond_num_est = 1;
|
|
unsigned freq_of_CN_estimation = 1;
|
|
|
|
Vec<dbl> prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(prediction_result,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_cond_num_est,
|
|
freq_of_CN_estimation,
|
|
tracking_tolerance);
|
|
|
|
BOOST_CHECK(success_code == bertini::SuccessCode::MatrixSolveFailureFirstPartOfPrediction);
|
|
|
|
|
|
}
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(euler_predict_linear_algebra_fails_mp)
|
|
{
|
|
// Circle line homotopy has singular point at (x,y) = (1,-4) and t = .75
|
|
|
|
bertini::DefaultPrecision(TRACKING_TEST_MPFR_DEFAULT_DIGITS);
|
|
|
|
// Starting point in spacetime step
|
|
Vec<mpfr> current_space(2);
|
|
current_space << mpfr("1.0"), mpfr("-4.0");
|
|
|
|
// Starting time
|
|
mpfr current_time(".7500000000000000000000000001");
|
|
// Time step
|
|
mpfr delta_t("-0.1");
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x - 5*y) );
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
AMP.coefficient_bound = 5;
|
|
|
|
double tracking_tolerance = 1e-5;
|
|
double condition_number_estimate;
|
|
|
|
unsigned num_steps_since_last_cond_num_est = 1;
|
|
unsigned freq_of_CN_estimation = 1;
|
|
|
|
Vec<mpfr> prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(prediction_result,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_cond_num_est,
|
|
freq_of_CN_estimation,
|
|
tracking_tolerance);
|
|
|
|
BOOST_CHECK(success_code == bertini::SuccessCode::MatrixSolveFailureFirstPartOfPrediction);
|
|
}
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(euler_predict_linear_criterion_a_is_false_d)
|
|
{
|
|
// Circle line homotopy has singular point at (x,y) = (1,-4) and t = .75
|
|
|
|
// Starting point in spacetime step
|
|
Vec<dbl> current_space(2);
|
|
current_space << dbl(1.0), dbl(-4.0);
|
|
|
|
// Starting time
|
|
dbl current_time(.8);
|
|
// Time step
|
|
dbl delta_t(-0.1);
|
|
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x - 5*y) );
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
double norm_J, norm_J_inverse, size_proportion, error_est;
|
|
|
|
AMP.coefficient_bound = 5;
|
|
AMP.safety_digits_1 = 100;
|
|
|
|
double tracking_tolerance(1e-5);
|
|
double condition_number_estimate;
|
|
|
|
unsigned num_steps_since_last_cond_num_est = 1;
|
|
unsigned freq_of_CN_estimation = 1;
|
|
|
|
Vec<dbl> prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(prediction_result,
|
|
size_proportion,
|
|
norm_J, norm_J_inverse,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_cond_num_est,
|
|
freq_of_CN_estimation,
|
|
tracking_tolerance,
|
|
AMP);
|
|
|
|
BOOST_CHECK(success_code == bertini::SuccessCode::HigherPrecisionNecessary);
|
|
}
|
|
|
|
BOOST_AUTO_TEST_CASE(euler_predict_linear_criterion_a_is_false_mp)
|
|
{
|
|
// Circle line homotopy has singular point at (x,y) = (1,-4) and t = .75
|
|
bertini::DefaultPrecision(TRACKING_TEST_MPFR_DEFAULT_DIGITS);
|
|
|
|
// Starting point in spacetime step
|
|
Vec<mpfr> current_space(2);
|
|
current_space << mpfr("1.0"), mpfr("-4.0");
|
|
|
|
// Starting time
|
|
mpfr current_time(".8");
|
|
// Time step
|
|
mpfr delta_t("-0.1");
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x - 5*y) );
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
double norm_J, norm_J_inverse, size_proportion, error_est;
|
|
|
|
AMP.coefficient_bound = 5;
|
|
AMP.safety_digits_1 = 100;
|
|
|
|
double tracking_tolerance = 1e-5;
|
|
double condition_number_estimate;
|
|
|
|
unsigned num_steps_since_last_cond_num_est = 1;
|
|
unsigned freq_of_CN_estimation = 1;
|
|
|
|
Vec<mpfr> prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(prediction_result,
|
|
size_proportion,
|
|
norm_J, norm_J_inverse,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_cond_num_est,
|
|
freq_of_CN_estimation,
|
|
tracking_tolerance,
|
|
AMP);
|
|
|
|
BOOST_CHECK(success_code == bertini::SuccessCode::HigherPrecisionNecessary);
|
|
}
|
|
|
|
BOOST_AUTO_TEST_CASE(euler_predict_linear_criterion_c_is_false_d)
|
|
{
|
|
// Circle line homotopy has singular point at (x,y) = (1,-4) and t = .75
|
|
|
|
// Starting point in spacetime step
|
|
Vec<dbl> current_space(2);
|
|
current_space << dbl(1.0), dbl(-4.0);
|
|
|
|
// Starting time
|
|
dbl current_time(.8);
|
|
// Time step
|
|
dbl delta_t(-0.1);
|
|
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x - 5*y) );
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
double norm_J, norm_J_inverse, size_proportion, error_est;
|
|
|
|
AMP.coefficient_bound = 5;
|
|
AMP.safety_digits_2 = 100;
|
|
|
|
AMP.SetPhiPsiFromBounds();
|
|
|
|
double tracking_tolerance(1e-5);
|
|
double condition_number_estimate;
|
|
|
|
unsigned num_steps_since_last_cond_num_est = 1;
|
|
unsigned freq_of_CN_estimation = 1;
|
|
|
|
Vec<dbl> prediction_result;
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(prediction_result,
|
|
size_proportion,
|
|
norm_J, norm_J_inverse,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_cond_num_est,
|
|
freq_of_CN_estimation,
|
|
tracking_tolerance,
|
|
AMP);
|
|
|
|
BOOST_CHECK(success_code == bertini::SuccessCode::HigherPrecisionNecessary);
|
|
}
|
|
|
|
BOOST_AUTO_TEST_CASE(euler_predict_linear_criterion_c_is_false_mp)
|
|
{
|
|
// Circle line homotopy has singular point at (x,y) = (1,-4) and t = .75
|
|
|
|
bertini::DefaultPrecision(TRACKING_TEST_MPFR_DEFAULT_DIGITS);
|
|
// Starting point in spacetime step
|
|
Vec<mpfr> current_space(2);
|
|
current_space << mpfr("1.0"), mpfr("-4.0");
|
|
|
|
// Starting time
|
|
mpfr current_time(".8");
|
|
// Time step
|
|
mpfr delta_t("-0.1");
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x - 5*y) );
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
double norm_J, norm_J_inverse, size_proportion, error_est;
|
|
|
|
AMP.coefficient_bound = 5;
|
|
AMP.safety_digits_2 = 100;
|
|
|
|
double tracking_tolerance = 1e-5;
|
|
double condition_number_estimate;
|
|
|
|
unsigned num_steps_since_last_cond_num_est = 1;
|
|
unsigned freq_of_CN_estimation = 1;
|
|
|
|
Vec<mpfr> prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(prediction_result,
|
|
size_proportion,
|
|
norm_J, norm_J_inverse,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_cond_num_est,
|
|
freq_of_CN_estimation,
|
|
tracking_tolerance,
|
|
AMP);
|
|
|
|
BOOST_CHECK(success_code == bertini::SuccessCode::HigherPrecisionNecessary);
|
|
}
|
|
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(circle_line_euler_change_precision)
|
|
{
|
|
bertini::DefaultPrecision(TRACKING_TEST_MPFR_DEFAULT_DIGITS);
|
|
|
|
// Starting point in spacetime step
|
|
Vec<mpfr> current_space(2);
|
|
current_space << mpfr("2.3","0.2"), mpfr("1.1", "1.87");
|
|
|
|
// Starting time
|
|
mpfr current_time("0.9");
|
|
// Time step
|
|
mpfr delta_t("-0.1");
|
|
|
|
|
|
|
|
|
|
bertini::System sys;
|
|
Var x = Variable::Make("x"), y = Variable::Make("y"), t = Variable::Make("t");
|
|
|
|
VariableGroup vars{x,y};
|
|
|
|
sys.AddVariableGroup(vars);
|
|
sys.AddPathVariable(t);
|
|
|
|
// Define homotopy system
|
|
sys.AddFunction( t*(pow(x,2)-1) + (1-t)*(pow(x,2) + pow(y,2) - 4) );
|
|
sys.AddFunction( t*(y-1) + (1-t)*(2*x + 5*y) );
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
BOOST_CHECK_EQUAL(AMP.degree_bound,2);
|
|
AMP.coefficient_bound = 5;
|
|
|
|
double norm_J, norm_J_inverse, size_proportion;
|
|
|
|
Vec<mpfr> predicted(2);
|
|
predicted << mpfr("2.40310963516214640018253210912048","0.187706567388887830930493342816564"),
|
|
mpfr("0.370984337833979085688209698697074", "1.30889906180158745272421523674049");
|
|
|
|
Vec<mpfr> euler_prediction_result;
|
|
|
|
double tracking_tolerance = 1e-5;
|
|
double condition_number_estimate;
|
|
unsigned num_steps_since_last_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::Euler,sys);
|
|
|
|
auto success_code = predictor->Predict(euler_prediction_result,
|
|
size_proportion,
|
|
norm_J, norm_J_inverse,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_condition_number_computation,
|
|
frequency_of_CN_estimation,
|
|
tracking_tolerance,
|
|
AMP);
|
|
|
|
BOOST_CHECK(success_code==bertini::SuccessCode::Success);
|
|
BOOST_CHECK_EQUAL(euler_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < euler_prediction_result.size(); ++ii)
|
|
BOOST_CHECK(abs(euler_prediction_result(ii)-predicted(ii)) < threshold_clearance_mp);
|
|
|
|
|
|
|
|
|
|
bertini::DefaultPrecision(50);
|
|
Precision(current_space,50); assert(current_space(0).precision()==50);
|
|
current_time.precision(50);
|
|
delta_t.precision(50);
|
|
sys.precision(50);assert(sys.precision()==50);
|
|
Precision(predicted,50); assert(predicted(0).precision()==50);
|
|
Precision(euler_prediction_result,50); assert(euler_prediction_result(0).precision()==50);
|
|
predictor->ChangePrecision(50);
|
|
|
|
// Starting point in spacetime step
|
|
current_space << mpfr("2.3","0.2"), mpfr("1.1", "1.87");
|
|
|
|
// Starting time
|
|
current_time = mpfr("0.9");
|
|
// Time step
|
|
delta_t = mpfr("-0.1");
|
|
|
|
|
|
|
|
|
|
|
|
|
|
AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
BOOST_CHECK_EQUAL(AMP.degree_bound,2);
|
|
AMP.coefficient_bound = 5;
|
|
|
|
|
|
predicted << mpfr("2.4031096351621464001825321091204810500227008230702","0.18770656738888783093049334281656388282191769240875"),
|
|
mpfr("0.37098433783397908568820969869707413571104273956141", "1.3088990618015874527242152367404908738825831867988");
|
|
|
|
|
|
success_code = predictor->Predict(euler_prediction_result,
|
|
size_proportion,
|
|
norm_J, norm_J_inverse,
|
|
sys,
|
|
current_space, current_time,
|
|
delta_t,
|
|
condition_number_estimate,
|
|
num_steps_since_last_condition_number_computation,
|
|
frequency_of_CN_estimation,
|
|
tracking_tolerance,
|
|
AMP);
|
|
|
|
BOOST_CHECK(success_code==bertini::SuccessCode::Success);
|
|
BOOST_CHECK_EQUAL(euler_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < euler_prediction_result.size(); ++ii)
|
|
{
|
|
BOOST_CHECK(abs(euler_prediction_result(ii)-predicted(ii)) < 1e-47);
|
|
}
|
|
}
|
|
|
|
BOOST_AUTO_TEST_SUITE_END()
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|