added first version of triple doubles

libmandel/include/Generators.h

@@ -22,6 +22,7 @@ namespace mnd
         DOUBLE_FLOAT,
         DOUBLE,
         DOUBLE_DOUBLE,
+        TRIPLE_DOUBLE,
         FLOAT128,
         FLOAT256,
         FLOAT512,
@@ -55,6 +56,7 @@ namespace mnd
     template<> Real getPrecision<float>();
     template<> Real getPrecision<double>();
     template<> Real getPrecision<DoubleDouble>();
+    template<> Real getPrecision<TripleDouble>();
     template<> Real getPrecision<QuadDouble>();
     template<> Real getPrecision<Fixed64>();
     template<> Real getPrecision<Fixed128>();
@@ -68,6 +70,7 @@ namespace mnd
     template<> inline Precision getType<float>() { return Precision::FLOAT; }
     template<> inline Precision getType<double>() { return Precision::DOUBLE; }
     template<> inline Precision getType<DoubleDouble>() { return Precision::DOUBLE_DOUBLE; }
+    template<> inline Precision getType<TripleDouble>() { return Precision::TRIPLE_DOUBLE; }
     template<> inline Precision getType<QuadDouble>() { return Precision::QUAD_DOUBLE; }
     template<> inline Precision getType<Fixed64>() { return Precision::FIXED64; }
     template<> inline Precision getType<Fixed128>() { return Precision::FIXED128; }
@@ -100,6 +103,7 @@ enum class mnd::GeneratorType : int
     DOUBLE_DOUBLE_AVX,
     DOUBLE_DOUBLE_AVX_FMA,
     DOUBLE_DOUBLE_NEON,
+    TRIPLE_DOUBLE,
     QUAD_DOUBLE,
     QUAD_DOUBLE_AVX_FMA,
     FLOAT128,

libmandel/include/PolyfloatUtil.h

@@ -0,0 +1,74 @@
+#ifndef MANDEL_POLYFLOATUTIL_H
+#define MANDEL_POLYFLOATUTIL_H
+
+namespace mnd
+{
+    namespace pfu
+    {
+        template<typename T>
+        struct Pair
+        {
+            T a;
+            T b;
+        };
+
+        using DoublePair = Pair<double>;
+        using FloatPair = Pair<float>;
+
+        template<typename T>
+        inline static Pair<T> quickTwoSum(T a, T b)
+        {
+            T s = a + b;
+            T e = b - (s - a);
+            return { s, e };
+        }
+
+        template<typename T>
+        inline Pair<T> twoSum(T a, T b)
+        {
+            T s = a + b;
+            T v = s - a;
+            T e = (a - (s - v)) + (b - v);
+            return { s, e };
+        }
+
+        template<typename T>
+        inline Pair<T> threeTwoSum(T a, T b, T c)
+        {
+            auto[t1, t2] = twoSum(a, b);
+            auto[r0, t3] = twoSum(t1, c);
+            return { r0, t2 + t3 };
+        }
+
+        inline Pair<double> split(double a)
+        {
+            static const double splitter = double((1ULL << 27) + 1);
+            double t = splitter * a;
+            double ahi = t - (t - a);
+            double alo = a - ahi;
+            return { ahi, alo };
+        }
+
+        template<typename T>
+        inline Pair<T> twoProd(T a, T b)
+        {
+            T p = a * b;
+            auto [ahi, alo] = split(a);
+            auto [bhi, blo] = split(b);
+            T e = ((ahi * bhi - p) + ahi * blo + alo * bhi) + alo * blo;
+            return { p, e };
+        }
+
+        template<typename T>
+        inline Pair<T> twoProdFma(T a, T b)
+        {
+            T p = a * b;
+            //T e = std::fma(a, b, -p);
+            T e = 0;
+            return { p, e };
+        }
+    }
+}
+
+#endif // MANDEL_POLYFLOATUTIL_H
+

libmandel/include/TripleDouble.h

@@ -0,0 +1,85 @@
+#ifndef MANDEL_TRIPLEDOUBLE_H
+#define MANDEL_TRIPLEDOUBLE_H
+
+#include "PolyfloatUtil.h"
+
+namespace mnd
+{
+    struct TripleDouble;
+}
+
+struct mnd::TripleDouble
+{
+    double x[3];
+
+    inline TripleDouble(double val) :
+        x{ val, 0.0, 0.0 }
+    {}
+
+    inline TripleDouble(double a, double b, double c) :
+        x{ a, b, c }
+    {}
+
+    double operator[]           (int i) const   { return x[i]; }
+    const double& operator[]    (int i)         { return x[i]; }
+};
+
+
+inline mnd::TripleDouble operator+(const mnd::TripleDouble& a,
+    const mnd::TripleDouble& b)
+{
+    auto[r0, t0] = mnd::pfu::twoSum(a[0], b[0]);
+    auto[t1, t2] = mnd::pfu::twoSum(a[1], b[1]);
+    auto[r1, t3] = mnd::pfu::twoSum(t0, t1);
+    auto r2 = t2 + t3 + a[2] + b[2];
+
+    auto[re1, t4] = mnd::pfu::quickTwoSum(r0, r1);
+    auto[re2, re3] = mnd::pfu::quickTwoSum(t4, r2);
+    return { re1, re2, re3 };
+}
+
+
+inline bool operator>(const mnd::TripleDouble& a, const mnd::TripleDouble& b)
+{
+    if (a[0] > b[0])
+        return true;
+    if (a[0] < b[0])
+        return false;
+    if (a[1] > b[1])
+        return true;
+    if (a[1] < b[1])
+        return false;
+    return a[2] > b[2];
+}
+
+inline mnd::TripleDouble operator-(const mnd::TripleDouble& a,
+    const mnd::TripleDouble& b)
+{
+    auto[r0, t0] = mnd::pfu::twoSum(a[0], -b[0]);
+    auto[t1, t2] = mnd::pfu::twoSum(a[1], -b[1]);
+    auto[r1, t3] = mnd::pfu::twoSum(t0, t1);
+    auto r2 = t2 + t3 + a[2] - b[2];
+
+    auto[re1, t4] = mnd::pfu::quickTwoSum(r0, r1);
+    auto[re2, re3] = mnd::pfu::quickTwoSum(t4, r2);
+    return { re1, re2, re3 };
+}
+
+
+inline mnd::TripleDouble operator*(const mnd::TripleDouble& a,
+    const mnd::TripleDouble& b)
+{
+    auto[p1_0, p2_0] = mnd::pfu::twoProd(a[0], b[0]);
+    auto[p2_1, p3_0] = mnd::pfu::twoProd(a[0], b[1]);
+    auto[p2_2, p3_1] = mnd::pfu::twoProd(a[1], b[0]);
+
+    auto[t2, tl3] = mnd::pfu::threeTwoSum(p2_0, p2_1, p2_2);
+    auto t3 = tl3 + p3_0 + p3_1 + a[1] * b[1] + a[2] * b[0] + a[0] * b[2];
+    auto[re0, q2] = mnd::pfu::quickTwoSum(p1_0, t2);
+    auto[re1, re2] = mnd::pfu::quickTwoSum(q2, t3);
+    return { re0, re1, re2 };
+}
+
+
+#endif // MANDEL_TRIPLEDOUBLE_H
+

libmandel/include/Types.h

@@ -23,6 +23,7 @@
 #include <qd/qd_real.h>
 
 #include "LightDoubleDouble.h"
+#include "TripleDouble.h"
 
 namespace mnd
 {
@@ -190,6 +191,35 @@ namespace mnd
     }
 
     template<>
+    inline Real convert<Real, LightDoubleDouble>(const LightDoubleDouble& x)
+    {
+        return Real{ x[0] } + x[1];
+    }
+
+    template<>
+    inline TripleDouble convert<TripleDouble, Real>(const Real& x)
+    {
+        double s = static_cast<double>(x);
+        Real tmp = x - s;
+        double e1 = static_cast<double>(tmp);
+        tmp = tmp - e1;
+        double e2 = static_cast<double>(tmp);
+        return TripleDouble { s, e1, e2 };
+    }
+
+    template<>
+    inline Real convert<Real, TripleDouble>(const TripleDouble& x)
+    {
+        return Real{ x[0] } + Real{ x[1] } + Real{ x[2] };
+    }
+
+    template<>
+    inline float convert<float, TripleDouble>(const TripleDouble& x)
+    {
+        return float(x.x[0] + x.x[1] + x.x[2]);
+    }
+
+    template<>
     inline QuadDouble convert<QuadDouble, Real>(const Real& x)
     {
         double s = static_cast<double>(x);

libmandel/src/CpuGenerators.cpp

@@ -48,6 +48,9 @@ namespace mnd
     template class CpuGenerator<mnd::DoubleDouble, mnd::NONE, false>;
     template class CpuGenerator<mnd::DoubleDouble, mnd::NONE, true>;
 
+    template class CpuGenerator<mnd::TripleDouble, mnd::NONE, false>;
+    template class CpuGenerator<mnd::TripleDouble, mnd::NONE, true>;
+
     template class CpuGenerator<mnd::QuadDouble, mnd::NONE, false>;
     template class CpuGenerator<mnd::QuadDouble, mnd::NONE, true>;
 #endif

libmandel/src/Generators.cpp

@@ -172,6 +172,7 @@ namespace mnd
             { Precision::DOUBLE_FLOAT, Real("4.0e-15") },
             { Precision::DOUBLE, getPrecision<double>() },
             { Precision::DOUBLE_DOUBLE, Real("1.0e-29") },
+            { Precision::TRIPLE_DOUBLE, Real("1.0e-45") },
             { Precision::QUAD_DOUBLE, Real("1.0e-56") },
             { Precision::FIXED64, Real("3.5e-15") },
             { Precision::FIXED128, Real("1.317e-29") },
@@ -208,6 +209,8 @@ namespace mnd
         case GeneratorType::DOUBLE_DOUBLE_AVX:
         case GeneratorType::DOUBLE_DOUBLE_AVX_FMA:
             return getPrecision<DoubleDouble>();
+        case GeneratorType::TRIPLE_DOUBLE:
+            return getPrecision<TripleDouble>();
         case GeneratorType::QUAD_DOUBLE:
             return getPrecision<QuadDouble>();
         case GeneratorType::FLOAT128:
@@ -241,6 +244,10 @@ namespace mnd
         return Real("1.0e-29");
     }
     template<>
+    Real getPrecision<TripleDouble>() {
+        return Real("1.0e-46");
+    }
+    template<>
     Real getPrecision<QuadDouble>() {
         return Real("3.0e-64");
     }

libmandel/src/Mandel.cpp

@@ -47,6 +47,7 @@ static const std::map<mnd::GeneratorType, std::string> typeNames =
     { mnd::GeneratorType::DOUBLE_DOUBLE_AVX, "double double AVX" },
     { mnd::GeneratorType::DOUBLE_DOUBLE_AVX_FMA, "double double AVX+FMA" },
     { mnd::GeneratorType::DOUBLE_DOUBLE_NEON, "double double NEON" },
+    { mnd::GeneratorType::TRIPLE_DOUBLE, "triple double" },
     { mnd::GeneratorType::QUAD_DOUBLE, "quad double" },
     { mnd::GeneratorType::QUAD_DOUBLE_AVX_FMA, "quad double AVX+FMA" },
     { mnd::GeneratorType::FLOAT128, "float128" },
@@ -192,6 +193,9 @@ MandelContext::MandelContext(void)
     cpuGenerators.insert({ GeneratorType::DOUBLE_DOUBLE, std::move(dd) });
     cpuGenerators.insert({ GeneratorType::QUAD_DOUBLE, std::move(qd) });
 
+    auto td = std::make_unique<CpuGenerator<TripleDouble, mnd::NONE, true>>();
+    cpuGenerators.insert({ GeneratorType::TRIPLE_DOUBLE, std::move(td) });
+
     auto fix512 = std::make_unique<CpuGenerator<Fixed512, mnd::NONE, true>>();
     cpuGenerators.insert({ GeneratorType::FIXED512, std::move(fix512) });
 

resources/Almond.qrc

@@ -13,4 +13,7 @@
         <file alias="peach">gradients/peach.xml</file>
         <file alias="blue gold">gradients/blue_gold.xml</file>
     </qresource>
+    <qresource prefix="/about">
+        <file alias="about_en">about_en.html</file>
+    </qresource>
 </RCC>