790 lines
22 KiB
C++
790 lines
22 KiB
C++
//This file is part of Bertini 2.0.
|
|
//
|
|
//heun_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.
|
|
//
|
|
//heun_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 heun_test.cpp. If not, see <http://www.gnu.org/licenses/>.
|
|
//
|
|
|
|
// Heun_test.cpp
|
|
//
|
|
// copyright 2016
|
|
// James B. Collins
|
|
// West Texas A&M University
|
|
// Department of Mathematics
|
|
// Spring 2016
|
|
|
|
|
|
|
|
#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(heun_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>;
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(circle_line_heun_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, error_est;
|
|
|
|
Vec<dbl> predicted(2);
|
|
predicted << dbl(2.38948874619536140814029774733947,0.208678935223681033727262214382917),
|
|
dbl(0.524558056401030798191044945035673, 1.43029356995029310361616395235936);
|
|
double predicted_error = .197349645229023708608160063982175;
|
|
|
|
Vec<dbl> heun_prediction_result;
|
|
dbl next_time;
|
|
|
|
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::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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(heun_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < heun_prediction_result.size(); ++ii)
|
|
{
|
|
BOOST_CHECK(abs(heun_prediction_result(ii)-predicted(ii)) < threshold_clearance_d);
|
|
}
|
|
BOOST_CHECK(fabs(error_est / predicted_error - 1) < threshold_clearance_d);
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(circle_line_heun_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, error_est;
|
|
|
|
Vec<mpfr> predicted(2);
|
|
predicted << mpfr("2.38948874619536140814029774733947","0.208678935223681033727262214382917"),
|
|
mpfr("0.524558056401030798191044945035673", "1.43029356995029310361616395235936");
|
|
double predicted_error = double(.197349645229023708608160063982175);
|
|
|
|
Vec<mpfr> heun_prediction_result;
|
|
mpfr next_time;
|
|
|
|
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::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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(heun_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < heun_prediction_result.size(); ++ii)
|
|
BOOST_CHECK(abs(heun_prediction_result(ii)-predicted(ii)) < threshold_clearance_mp);
|
|
|
|
using std::abs;
|
|
BOOST_CHECK(abs(error_est - predicted_error) < std::numeric_limits<double>::epsilon());
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(monodromy_heun_d)
|
|
{
|
|
bertini::DefaultPrecision(TRACKING_TEST_MPFR_DEFAULT_DIGITS);
|
|
|
|
// Starting point in spacetime step
|
|
Vec<dbl> current_space(2);
|
|
current_space << dbl(0.464158883361277585510862309093), dbl(0.74161984870956629487113974408);
|
|
|
|
// 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");
|
|
std::shared_ptr<Float> half = Float::Make("0.5");
|
|
|
|
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) + half) );
|
|
|
|
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
double norm_J, norm_J_inverse, size_proportion, error_est;
|
|
|
|
BOOST_CHECK_EQUAL(AMP.degree_bound,3);
|
|
AMP.coefficient_bound = 2;
|
|
|
|
|
|
Vec<dbl> predicted(2);
|
|
predicted << dbl(0.412299156269677938503694812160886),
|
|
dbl(0.731436945256924470273568899877140);
|
|
double predicted_error = 0.00544428757292458409463632380167773;
|
|
|
|
Vec<dbl> heun_prediction_result;
|
|
double next_time;
|
|
|
|
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::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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(heun_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < heun_prediction_result.size(); ++ii)
|
|
{
|
|
BOOST_CHECK(abs(heun_prediction_result(ii)-predicted(ii)) < threshold_clearance_d);
|
|
}
|
|
|
|
BOOST_CHECK(fabs(error_est / predicted_error - 1) < threshold_clearance_d);
|
|
|
|
}
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(monodromy_heun_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");
|
|
std::shared_ptr<Float> half = Float::Make("0.5");
|
|
|
|
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) + half) );
|
|
|
|
|
|
|
|
auto AMP = bertini::tracking::AMPConfigFrom(sys);
|
|
|
|
BOOST_CHECK_EQUAL(AMP.degree_bound,3);
|
|
AMP.coefficient_bound = 2;
|
|
|
|
double norm_J, norm_J_inverse, size_proportion, error_est;
|
|
|
|
|
|
Vec<mpfr> predicted(2);
|
|
predicted << mpfr("0.412299156269677938503694812160886"),
|
|
mpfr("0.731436945256924470273568899877140");
|
|
double predicted_error = double(0.00544428757292458409463632380167773);
|
|
|
|
Vec<mpfr> heun_prediction_result;
|
|
mpfr next_time;
|
|
|
|
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::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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(heun_prediction_result.size(),2);
|
|
for (unsigned ii = 0; ii < heun_prediction_result.size(); ++ii)
|
|
{
|
|
BOOST_CHECK(abs(heun_prediction_result(ii)-predicted(ii)) < threshold_clearance_mp);
|
|
}
|
|
|
|
using std::abs;
|
|
BOOST_CHECK(abs(error_est / predicted_error - 1) < std::numeric_limits<double>::epsilon());
|
|
}
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(heun_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);
|
|
// 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;
|
|
|
|
double tracking_tolerance(1e-5);
|
|
double condition_number_estimate;
|
|
|
|
unsigned num_steps_since_last_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
Vec<dbl> heun_prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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::MatrixSolveFailureFirstPartOfPrediction);
|
|
|
|
}
|
|
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(heun_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);
|
|
|
|
double norm_J, norm_J_inverse, size_proportion, error_est;
|
|
|
|
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> heun_prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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::MatrixSolveFailureFirstPartOfPrediction);
|
|
}
|
|
|
|
|
|
BOOST_AUTO_TEST_CASE(heun_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_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
Vec<dbl> heun_prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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::HigherPrecisionNecessary);
|
|
}
|
|
|
|
BOOST_AUTO_TEST_CASE(heun_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_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
Vec<mpfr> heun_prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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::HigherPrecisionNecessary);
|
|
}
|
|
|
|
BOOST_AUTO_TEST_CASE(heun_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_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
Vec<dbl> heun_prediction_result;
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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::HigherPrecisionNecessary);
|
|
}
|
|
|
|
BOOST_AUTO_TEST_CASE(heun_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_condition_number_computation = 1;
|
|
unsigned frequency_of_CN_estimation = 1;
|
|
|
|
Vec<mpfr> heun_prediction_result;
|
|
|
|
|
|
std::shared_ptr<ExplicitRKPredictor> predictor = std::make_shared< ExplicitRKPredictor >(bertini::tracking::Predictor::HeunEuler,sys);
|
|
|
|
auto success_code = predictor->Predict(heun_prediction_result,
|
|
error_est,
|
|
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::HigherPrecisionNecessary);
|
|
}
|
|
|
|
|
|
BOOST_AUTO_TEST_SUITE_END()
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|