#include <algorithm>
#include <cmath>
#include <iostream>
#include <limits>
#include <stdexcept>
bool close(double value,double reference,double atol=1e-12,double rtol=1e-9) {
    if(!std::isfinite(value)||!std::isfinite(reference)||!std::isfinite(atol)||!std::isfinite(rtol)||atol<0||rtol<0) return false;
    const double difference=std::abs(value-reference);
    const double tolerance=atol+rtol*std::abs(reference);
    if(std::isfinite(difference)) return difference<=tolerance;
    if(std::isfinite(tolerance)) return false;
    // Both overflowed: compare the same quantities in bounded units.
    // Do not rely on long double having extra range on this platform.
    const double scale=std::max(std::abs(value),std::abs(reference));
    return std::abs(value/scale-reference/scale)
        <= atol/scale+rtol*(std::abs(reference)/scale);
}
int main() {
    const double large=1e16,negative=-1e16,small=1;
    const double serial=(large+negative)+small,regrouped=large+(negative+small);
    const double nan=std::numeric_limits<double>::quiet_NaN(),inf=std::numeric_limits<double>::infinity();
    if(serial!=1 || regrouped!=0 || close(nan,1) || close(1,nan) || close(inf,inf)
       || !close(1+1e-10,1) || !close(1e-13,0) || close(1,0) || close(1,1,-1))
        throw std::runtime_error("numeric contract");
    const double maximum=std::numeric_limits<double>::max();
    if(close(maximum,-maximum/2,0,2.5) || !close(maximum,-maximum/2,0,3)
       || close(maximum,-maximum,0,1.5) || !close(maximum,-maximum,0,2)
       || close(maximum,-maximum,maximum,0) || !close(maximum,-maximum,maximum,1)
       || !close(maximum,maximum/2,0,2) || !close(maximum,maximum,0,0))
        throw std::runtime_error("finite extreme tolerance");
    std::cout<<"serial="<<serial<<" regrouped="<<regrouped<<" nonfinite_rejected=1\n";
}
