fix(endgames): NaN is a failure, never Converged; NaN-aware security valve

core/include/bertini2/eigen_extensions.hpp

@@ -312,6 +312,26 @@ namespace bertini {
 	}
 
 
+	/**
+	\brief True when any entry of a complex-valued Eigen object has a NaN real or imaginary part.
+
+	Component-wise isnan is REQUIRED for correctness: the self-inequality trick (`z != z`) does
+	not work for complex_mp, whose equality does not follow IEEE NaN semantics (a NaN complex_mp
+	compares EQUAL to itself) -- which also makes Eigen's own hasNaN() unreliable there.
+	*/
+	template<typename Derived>
+	inline
+	bool ContainsNaN(Eigen::MatrixBase<Derived> const & v)
+	{
+		using std::isnan;
+		for (Eigen::Index ii = 0; ii < v.rows(); ++ii)
+			for (Eigen::Index jj = 0; jj < v.cols(); ++jj)
+				if (isnan(v(ii,jj).real()) || isnan(v(ii,jj).imag()))
+					return true;
+		return false;
+	}
+
+
 	/**
 	\brief Get the precision of an Eigen object. If the object is empty, it's the precision of a default-constructed Scalar.  If it actually has content, then it's the precision of the first element.
 

core/include/bertini2/endgames/base_endgame.hpp

@@ -326,6 +326,17 @@ class EndgameBase :
 		return this->template Get<SecurityConfig>();
 	}
 
+	/// \return True when a dehomogenized norm calls for security truncation: above max_norm, or
+	/// NaN.  NaN must count as beyond -- every IEEE comparison against NaN is false, so a plain
+	/// `norm > max_norm` valve is BLIND to exactly the paths most likely to be diverging (a NaN
+	/// dehom norm means the homogenizing coordinate vanished: the path is at infinity).
+	template<typename RealT>
+	bool BeyondSecurityMaxNorm(RealT const& norm) const
+	{
+		using std::isnan;
+		return isnan(norm) || norm > static_cast<RealT>(this->SecuritySettings().max_norm);
+	}
+
 	/// \brief Construct the endgame for a tracker, with its configuration as a tuple.
 	explicit EndgameBase(TrackerType const& tr, const ConfigsAsTuple& settings ) :
       	EndgamePrecPolicyBase<TrackerType>(tr), Configured( settings ), PrecT(tr)

core/include/bertini2/endgames/cauchy.hpp

@@ -1051,6 +1051,11 @@ class CauchyEndgame :
 			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;
 
 	}
@@ -1377,8 +1382,8 @@ class CauchyEndgame :
 				if (in_operating_zone)
 				{
 					norm_of_dehom_latest = this->GetSystem().InfinityNormOfDehomogenized(latest_approx);
-					if (norm_of_dehom_prev   > this->SecuritySettings().max_norm &&
-						norm_of_dehom_latest > this->SecuritySettings().max_norm  )
+					if (this->BeyondSecurityMaxNorm(norm_of_dehom_prev) &&
+						this->BeyondSecurityMaxNorm(norm_of_dehom_latest))
 					{
 						NotifyObservers(SecurityMaxNormReached<EmitterType>(*this));
 						return SuccessCode::SecurityMaxNormReached;
@@ -1515,8 +1520,8 @@ class CauchyEndgame :
 				if (in_operating_zone)
 				{
 					norm_of_dehom_latest = this->GetSystem().InfinityNormOfDehomogenized(this->final_approximation_);
-					if (norm_of_dehom_prev   > this->SecuritySettings().max_norm &&
-					    norm_of_dehom_latest > this->SecuritySettings().max_norm)
+					if (this->BeyondSecurityMaxNorm(norm_of_dehom_prev) &&
+					    this->BeyondSecurityMaxNorm(norm_of_dehom_latest))
 					{
 						NotifyObservers(SecurityMaxNormReached<EmitterType>(*this));
 						return SuccessCode::SecurityMaxNormReached;

core/include/bertini2/endgames/powerseries.hpp

@@ -443,6 +443,11 @@ class PowerSeriesEndgame :
 
 		auto min_found_difference = Eigen::NumTraits<RealT>::highest();
 
+		// Default the selection before the candidate loop: if every candidate's difference is
+		// NaN (poisoned samples), the comparisons below are all false and nothing is assigned --
+		// a fresh endgame would otherwise carry cycle number 0 into TransformToSPlane and throw.
+		this->cycle_number_ = 1;
+
 		TimeCont<ComplexT> s_times(num_pts);
 		SampCont<ComplexT> s_derivatives(num_pts);
 
@@ -589,6 +594,12 @@ class PowerSeriesEndgame :
 
 		Precision(result, Precision(s_derivatives.back()));
 		result = HermiteInterpolateAndSolve(ComplexT(0), num_pts, s_times, std::get<SampCont<ComplexT> >(samples_), s_derivatives, ContStart::Back);
+		// A NaN extrapolation must be a FAILURE code.  The run loop converges on
+		// `approx_error > FinalTolerance()` becoming false, and every IEEE comparison against
+		// NaN is false -- so a NaN approximation would exit the loop down the SUCCESS path,
+		// reporting Converged with a poisoned answer.
+		if (bertini::ContainsNaN(result))
+			return SuccessCode::FailedToConverge;
 		return SuccessCode::Success;
 	}//end ComputeApproximationOfXAtT0
 
@@ -783,7 +794,7 @@ class PowerSeriesEndgame :
 	 		if(this->SecuritySettings().level <= 0)
 	 		{
 	 			norm_of_dehom_of_latest_approx = this->GetSystem().InfinityNormOfDehomogenized(latest_approx);
-		 		if(norm_of_dehom_of_latest_approx > this->SecuritySettings().max_norm && norm_of_dehom_of_prev_approx > this->SecuritySettings().max_norm)
+		 		if(this->BeyondSecurityMaxNorm(norm_of_dehom_of_latest_approx) && this->BeyondSecurityMaxNorm(norm_of_dehom_of_prev_approx))
 		 		{
 		 			NotifyObservers(SecurityMaxNormReached<EmitterType>(*this));
 	 				return SuccessCode::SecurityMaxNormReached;
@@ -876,7 +887,7 @@ class PowerSeriesEndgame :
 			if (this->SecuritySettings().level <= 0)
 			{
 				norm_latest = this->GetSystem().InfinityNormOfDehomogenized(this->final_approximation_);
-				if (norm_latest > this->SecuritySettings().max_norm && norm_prev > this->SecuritySettings().max_norm)
+				if (this->BeyondSecurityMaxNorm(norm_latest) && this->BeyondSecurityMaxNorm(norm_prev))
 				{
 					NotifyObservers(SecurityMaxNormReached<EmitterType>(*this));
 					return SuccessCode::SecurityMaxNormReached;

core/test/endgames/generic_cauchy_test.hpp

@@ -1111,6 +1111,51 @@ BOOST_AUTO_TEST_CASE(compute_cauchy_approximation_cycle_num_1)
 }// end compute_cauchy_approximation_cycle_num_1
 
 
+/**
+Regression: a NaN in the loop samples must make the trapezoidal mean come back as a FAILURE
+code, never Success -- NaN compares false against everything, so downstream convergence and
+security comparisons are blind to a poisoned approximation.
+*/
+BOOST_AUTO_TEST_CASE(nan_cauchy_sample_yields_failure_code_not_success)
+{
+	DefaultPrecision(ambient_precision);
+
+	System sys;
+	Var x = Variable::Make("x");
+	Var t = Variable::Make("t");
+	sys.AddFunction((x-1)*(1-t) + (x+1)*t);
+	VariableGroup vars{x};
+	sys.AddVariableGroup(vars);
+	sys.AddPathVariable(t);
+
+	auto precision_config = PrecisionConfig(sys);
+	TrackerType tracker(sys);
+	bertini::tracking::SteppingConfig stepping_preferences;
+	bertini::tracking::NewtonConfig newton_preferences;
+	tracker.Setup(TestedPredictor, 1e-5, 1e5, stepping_preferences, newton_preferences);
+	tracker.PrecisionSetup(precision_config);
+
+	bertini::TimeCont<BCT> cauchy_times;
+	bertini::SampCont<BCT> cauchy_samples;
+	Vec<BCT> sample(1);
+
+	// cycle 1 x num_sample_points samples + the closing copy, one poisoned with NaN
+	cauchy_times.push_back(ComplexFromString("0.1"));        sample << ComplexFromString("0.8");  cauchy_samples.push_back(sample);
+	cauchy_times.push_back(ComplexFromString("0", "0.1"));   sample << ComplexFromString("0.81"); cauchy_samples.push_back(sample);
+	cauchy_times.push_back(ComplexFromString("-0.1"));       sample << BCT(std::numeric_limits<BRT>::quiet_NaN()); cauchy_samples.push_back(sample);
+	cauchy_times.push_back(ComplexFromString("0.1"));        sample << ComplexFromString("0.8");  cauchy_samples.push_back(sample);
+
+	TestedEGType my_endgame(tracker);
+	my_endgame.SetCauchyTimes(cauchy_times);
+	my_endgame.SetCauchySamples(cauchy_samples);
+	my_endgame.CycleNumber(1);
+
+	Vec<BCT> approx;
+	auto code = my_endgame.template ComputeCauchyApproximationOfXAtT0<BCT>(approx);
+	BOOST_CHECK(code != SuccessCode::Success);
+}// end nan_cauchy_sample_yields_failure_code_not_success
+
+
 
 /**
 	This test case uses all the sample points collected by CircleTrack around a non-singular point to compute an extrapolant using the 

core/test/endgames/generic_pseg_test.hpp

@@ -153,6 +153,55 @@ BOOST_AUTO_TEST_CASE( hermite_reproduces_low_degree_polynomials_exactly )
 }//end hermite_reproduces_low_degree_polynomials_exactly
 
 
+/**
+Regression: a NaN extrapolation must come back as a FAILURE code, never Success.  The run loop
+converges on `approx_error > FinalTolerance()` going false, and every IEEE comparison against NaN
+is false -- so a NaN approximation used to exit the loop down the SUCCESS path, reporting
+Converged with a poisoned answer.
+*/
+BOOST_AUTO_TEST_CASE(nan_sample_yields_failure_code_not_success)
+{
+	DefaultPrecision(ambient_precision);
+
+	bertini::System sys;
+	Var x = Variable::Make("x"), t = Variable::Make("t");
+	VariableGroup vars{x};
+	sys.AddVariableGroup(vars);
+	sys.AddPathVariable(t);
+	sys.AddFunction( pow(x-1,3)*(1-t) + (pow(x,3)+1)*t);
+
+	auto precision_config = PrecisionConfig(sys);
+	TrackerType tracker(sys);
+	bertini::tracking::SteppingConfig stepping_settings;
+	bertini::tracking::NewtonConfig newton_settings;
+	tracker.Setup(TestedPredictor, 1e-5, 1e5, stepping_settings, newton_settings);
+	tracker.PrecisionSetup(precision_config);
+
+	bertini::TimeCont<BCT> times;
+	bertini::SampCont<BCT> samples;
+	Vec<BCT> sample(1);
+
+	times.push_back(ComplexFromString(".1"));
+	sample << ComplexFromString("0.5");     samples.push_back(sample);
+	times.push_back(ComplexFromString(".05"));
+	sample << ComplexFromString("0.6");     samples.push_back(sample);
+	times.push_back(ComplexFromString(".025"));
+	sample << BCT(std::numeric_limits<BRT>::quiet_NaN());   samples.push_back(sample);   // poisoned
+
+	bertini::endgame::EndgameConfig endgame_settings;
+	TestedEGType my_endgame(tracker, endgame_settings);
+	my_endgame.SetTimes(times);
+	my_endgame.SetSamples(samples);
+	my_endgame.template ComputeAllDerivatives<BCT>();
+	my_endgame.SetRandVec<BCT>(1);
+	my_endgame.CycleNumber(1);
+
+	Vec<BCT> approx(1);
+	auto code = my_endgame.template ComputeApproximationOfXAtT0<BCT>(approx, BCT(0));
+	BOOST_CHECK(code != SuccessCode::Success);
+}//end nan_sample_yields_failure_code_not_success
+
+