Repair neg-zero and tiny-arg erf and erfc

[ckormanyos]

No description provided.

[jzmaddock]

Well done, I realized when I was lying in bed last night that 1 on it's own didn't quite cut it!

[ckormanyos]

It was a great review john. All you had to do was get me thinking. Thx.

[ckormanyos]

Actually John (@jzmaddock) in my new code, I think the final path for identically zero is not needed. That will be picked up by the first path, right?

[ckormanyos]

But then for the case of identically zoro, we would compute term, which is a waste of time.

Is this a case for (boost::math::tools::fpclassify)(z) since you check for NaN also?

[ckormanyos]

I am moving toward this implementation locally. I will also throw in a few more edge-case Multiprecision tests locally and wrap this one up shortly.

template <class T, class Policy, class Tag>
T erf_imp(T z, bool invert, const Policy& pol, const Tag& t)
{
   BOOST_MATH_STD_USING

   BOOST_MATH_INSTRUMENT_CODE("Generic erf_imp called");

   if ((boost::math::isnan)(z))
      return policies::raise_domain_error("boost::math::erf<%1%>(%1%)", "Expected a finite argument but got %1%", z, pol);

   const T abs_z { fabs(z) };

   const bool signbit_result = ((boost::math::signbit)(z) != 0);

   if (abs_z < tools::root_epsilon<T>())
   {
      if (abs_z < tools::epsilon<T>())
      {
         return invert ? T(1) : (!signbit_result) ? T(0) : -T(0);
      }

      // Series[Erf[x], {x, 0, 4}]
      // Series[Erfc[x], {x, 0, 4}]

      const T term2 { z * 2 / constants::root_pi<T>() };

      return invert ? 1 - term2 : term2;
   }

   if (signbit_result)
   {
      if(!invert)
         return -erf_imp(T(-z), invert, pol, t);
      else
         return 1 + erf_imp(T(-z), false, pol, t);
   }

   T result;

   if(!invert && (z > detail::erf_asymptotic_limit<T, Policy>()))
   {
      detail::erf_asympt_series_t<T> s(z);
      std::uintmax_t max_iter = policies::get_max_series_iterations<Policy>();
      result = boost::math::tools::sum_series(s, policies::get_epsilon<T, Policy>(), max_iter, 1);
      policies::check_series_iterations<T>("boost::math::erf<%1%>(%1%, %1%)", max_iter, pol);
   }
   else
   {
      T x = z * z;
      if(z < 1.3f)
      {
         // Compute P:
         // This is actually good for z p to 2 or so, but the cutoff given seems
         // to be the best compromise.  Performance wise, this is way quicker than anything else...
         result = erf_series_near_zero_sum(z, pol);
      }
      else if(x > 1 / tools::epsilon<T>())
      {
         // http://functions.wolfram.com/06.27.06.0006.02
         invert = !invert;
         result = exp(-x) / (constants::root_pi<T>() * z);
      }
      else
      {
         // Compute Q:
         invert = !invert;
         result = z * exp(-x);
         result /= boost::math::constants::root_pi<T>();
         result *= upper_gamma_fraction(T(0.5f), x, policies::get_epsilon<T, Policy>());
      }
   }

   return ((!invert) ? result : 1 - result);
}

[ckormanyos]

This is the generic implementation of erf/erfc in erf_imp. So it is used for MP types.

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 95.10%. Comparing base (4323902) to head (c1d2bec).
⚠️ Report is 9 commits behind head on develop.

Additional details and impacted files

@@           Coverage Diff            @@
##           develop    #1315   +/-   ##
========================================
  Coverage    95.10%   95.10%           
========================================
  Files          797      797           
  Lines        67118    67126    +8     
========================================
+ Hits         63833    63841    +8     
  Misses        3285     3285           

Files with missing lines	Coverage Δ
include/boost/math/special_functions/erf.hpp	100.00% <ø> (ø)
test/test_erf.hpp	100.00% <100.00%> (ø)

---

Continue to review full report in Codecov by Sentry.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 4323902...c1d2bec. Read the comment docs.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.