Nicolas Winkler 5 vuotta sitten
vanhempi
commit
ad59fc1977

+ 1 - 0
Almond.qrc

@@ -3,5 +3,6 @@
         <file alias="default">gradients/default.xml</file>
         <file alias="clouds">gradients/clouds.xml</file>
         <file alias="rainbow">gradients/rainbow.xml</file>
+        <file alias="grayscale">gradients/grayscale.xml</file>
     </qresource>
 </RCC>

+ 1 - 1
Gradient.cpp

@@ -86,7 +86,7 @@ Gradient Gradient::readXml(const QString& xml)
         uint8_t b = uint8_t(child.attributeNode("b").value().toInt());
         float p = child.attributeNode("p").value().toInt();
 
-        printf("rgb (%s): %d, %d, %d\n", child.text().toUtf8().data(), r, g, b);
+        //printf("rgb (%s): %d, %d, %d\n", child.text().toUtf8().data(), r, g, b);
         colorArr.push_back({ { r, g, b }, p });
     }
 

+ 10 - 15
MandelWidget.cpp

@@ -623,7 +623,7 @@ void MandelWidget::initializeGL(void)
     glDisable(GL_DEPTH_TEST);
 
     // looks not even better
-    //glDisable(GL_FRAMEBUFFER_SRGB);
+    glEnable(GL_FRAMEBUFFER_SRGB);
 
     //glShadeModel(GL_SMOOTH);
 
@@ -678,12 +678,13 @@ void MandelWidget::updateAnimations(void)
     else {
         auto now = std::chrono::high_resolution_clock::now();
         auto millis = std::chrono::duration_cast<std::chrono::milliseconds>(now - lastAnimUpdate).count();
-        double factor = ::pow(0.97, millis);
+        const mnd::Real factor = mnd::Real(::pow(0.97, millis));
+        const mnd::Real one(1.0);
 
-        currentViewport.x = currentViewport.x * factor + targetViewport.x * (1.0 - factor);
-        currentViewport.y = currentViewport.y * factor + targetViewport.y * (1.0 - factor);
-        currentViewport.width = currentViewport.width * factor + targetViewport.width * (1.0 - factor);
-        currentViewport.height = currentViewport.height * factor + targetViewport.height * (1.0 - factor);
+        currentViewport.x = currentViewport.x * factor + targetViewport.x * (one - factor);
+        currentViewport.y = currentViewport.y * factor + targetViewport.y * (one - factor);
+        currentViewport.width = currentViewport.width * factor + targetViewport.width * (one - factor);
+        currentViewport.height = currentViewport.height * factor + targetViewport.height * (one - factor);
 
         lastAnimUpdate = now;
         emit update();
@@ -720,7 +721,7 @@ void MandelWidget::drawInfo(void)
     const float DIST_FROM_BORDER = 15;
     float maxWidth = this->width() - 2 * DIST_FROM_BORDER;
     mnd::Real distPerPixel = currentViewport.width / this->width();
-    float log10 = (mnd::convert<float>(mnd::log(distPerPixel)) + ::logf(maxWidth)) / ::log(10);
+    float log10 = (mnd::convert<float>(mnd::log(distPerPixel)) + ::logf(maxWidth)) / ::logf(10);
     mnd::Real displayDist = mnd::pow(mnd::Real(10), ::floor(log10));
     float pixels = mnd::convert<float>(displayDist / distPerPixel);
     int factor = 1;
@@ -750,8 +751,8 @@ void MandelWidget::drawInfo(void)
     infoPainter.drawLine(QPointF{ DIST_FROM_BORDER, lineY }, QPointF{ lineXEnd, lineY });
     infoPainter.drawLine(QPointF{ DIST_FROM_BORDER, lineY }, QPointF{ DIST_FROM_BORDER, lineY - 5 });
     infoPainter.drawLine(QPointF{ lineXEnd, lineY }, QPointF{ lineXEnd, lineY - 5 });
-    infoPainter.drawText(DIST_FROM_BORDER, lineY - 20, lineXEnd - DIST_FROM_BORDER, 20, Qt::AlignCenter, QString::fromStdString(dis.str()));
-    //infoPainter.drawText(0, this->height() - 30, this->width(), 30, Qt::AlignBottom | Qt::AlignVCenter, QString::fromStdString(ss.str()));
+    infoPainter.drawText(int(DIST_FROM_BORDER), int(lineY - 20), int(lineXEnd - DIST_FROM_BORDER), 20,
+                         Qt::AlignCenter, QString::fromStdString(dis.str()));
     infoPainter.end();
 }
 
@@ -789,12 +790,6 @@ void MandelWidget::requestRecalc()
 }
 
 
-void MandelWidget::resizeGL(int width, int height)
-{
-    emit update();
-}
-
-
 void MandelWidget::resizeEvent(QResizeEvent* re)
 {
     QOpenGLWidget::resizeEvent(re);

+ 1 - 7
MandelWidget.h

@@ -116,6 +116,7 @@ struct PairHash {
     std::size_t operator () (const std::pair<T1, T2>& p) const {
         auto h1 = std::hash<T1>{}(p.first);
         auto h2 = std::hash<T2>{}(p.second);
+        //boost::hash_combine(h1, p.second);
         return (h1 ^ 234579245) * 23452354 + h2;
     }
 };
@@ -318,9 +319,6 @@ public:
     void setDisplayInfo(bool di);
 
     void initializeGL(void) override;
-
-    void resizeGL(int w, int h) override;
-
     void paintGL() override;
 
 private:
@@ -334,8 +332,6 @@ public:
     void setViewport(const mnd::MandelViewport& viewport);
     void setMaxIterations(int maxIter);
 
-    //void redraw();
-
     void requestRecalc(void);
 
     void resizeEvent(QResizeEvent* re) override;
@@ -347,7 +343,5 @@ public:
     inline const mnd::MandelViewport& getViewport(void) const { return targetViewport; }
 signals:
     void needsUpdate(const mnd::MandelInfo vp);
-public slots:
-    //void viewUpdated(Bitmap<RGBColor>* bitmap);
 };
 

+ 9 - 10
benchmarkdialog.cpp

@@ -5,7 +5,6 @@
 
 mnd::MandelViewport Benchmarker::benchViewport(void)
 {
-    //return mnd::MandelViewport{ -0.758267525104592591494, -0.066895616551111110830, 0.000000043217777777655, 0.000000043217777777655 };
     return mnd::MandelViewport{ -1.250000598933854152929, 0.0001879894057291665530, 0.0000003839916666666565, 0.0000003839916666666565 };
 }
 
@@ -60,7 +59,7 @@ std::pair<long long, std::chrono::nanoseconds> Benchmarker::measureMips(const st
 
     long long sum = 0;
     for (int i = 0; i < bitmap->width * bitmap->height; i++) {
-        sum += std::floor(bitmap->pixels[size_t(i)]);
+        sum += static_cast<long long>(std::floor(bitmap->pixels[size_t(i)]));
     }
 
     return std::make_pair(sum, duration_cast<nanoseconds>(after - before));
@@ -68,7 +67,7 @@ std::pair<long long, std::chrono::nanoseconds> Benchmarker::measureMips(const st
 
 double Benchmarker::benchmarkResult(mnd::Generator& mg) const
 {
-    int testIndex = 0;
+    size_t testIndex = 0;
 
     for (size_t i = 0; i < benches.size(); i++) {
         const mnd::MandelInfo& mi = benches[i];
@@ -205,14 +204,14 @@ BenchmarkDialog::BenchmarkDialog(mnd::MandelContext& mndContext, QWidget *parent
     printf("bench!\n"); fflush(stdout);
 
     auto& devices = mndContext.getDevices();
-    int nDevices = devices.size() + 1;
+    size_t nDevices = devices.size() + 1;
     ui.tableWidget->setColumnCount(6);
-    ui.tableWidget->setRowCount(nDevices);
+    ui.tableWidget->setRowCount(int(nDevices));
     ui.tableWidget->setHorizontalHeaderLabels({"Single Precision", "Double Precision", "Double-Double Precision", "Quad-Double Precision", "Quad Precision", "Oct Precision"});
 
     QString cpuDesc = ("CPU [" + mndContext.getCpuInfo().getBrand() + "]").c_str();
     ui.tableWidget->setVerticalHeaderItem(0, new QTableWidgetItem(cpuDesc));
-    for (int i = 0; i < devices.size(); i++) {
+    for (size_t i = 0; i < devices.size(); i++) {
         std::string cpuDescS = std::string("GPU ") + std::to_string(i + 1) + " [" + devices[i].getName().c_str() + "]";
         QString cpuDesc = QString::fromLatin1(cpuDescS.c_str());
         /*printf("brand [%d]: --> %s <--\n", (int) cpuDescS.size(), cpuDescS.c_str());
@@ -222,7 +221,7 @@ BenchmarkDialog::BenchmarkDialog(mnd::MandelContext& mndContext, QWidget *parent
         printf("\n");*/
         auto label = new QTableWidgetItem(cpuDesc);
         label->setStatusTip(QString::fromLatin1(devices[i].getName().c_str()));
-        ui.tableWidget->setVerticalHeaderItem(i + 1, label);
+        ui.tableWidget->setVerticalHeaderItem(int(i + 1), label);
     }
 
     qRegisterMetaType<BenchmarkResult>();
@@ -239,9 +238,9 @@ BenchmarkDialog::BenchmarkDialog(mnd::MandelContext& mndContext, QWidget *parent
 void BenchmarkDialog::update(BenchmarkResult br)
 {
     std::vector<double> cpu = br.values[0];
-    for (int j = 0; j < int(br.values.size()); j++) {
-        for (int i = 0; i < int(br.values[j].size()); i++) {
-            ui.tableWidget->setItem(j, i, new QTableWidgetItem(QString::number(br.values[j][i])));
+    for (size_t j = 0; j < br.values.size(); j++) {
+        for (size_t i = 0; i < br.values[j].size(); i++) {
+            ui.tableWidget->setItem(int(j), int(i), new QTableWidgetItem(QString::number(br.values[j][i])));
         }
     }
     ui.progressBar->setValue(int(br.percentage));

+ 2 - 2
benchmarks.ui

@@ -6,8 +6,8 @@
    <rect>
     <x>0</x>
     <y>0</y>
-    <width>864</width>
-    <height>401</height>
+    <width>1215</width>
+    <height>550</height>
    </rect>
   </property>
   <property name="windowTitle">

+ 1 - 0
gradientchoosedialog.cpp

@@ -15,6 +15,7 @@ GradientChooseDialog::GradientChooseDialog()
     gcd.presets->addItem("default");
     gcd.presets->addItem("clouds");
     gcd.presets->addItem("rainbow");
+    gcd.presets->addItem("grayscale");
 }
 
 

+ 1 - 2
libmandel/CMakeLists.txt

@@ -51,8 +51,7 @@ FILE(GLOB QdSources qd-2.3.22/src/*.cpp)
 
 target_compile_definitions(mandel PUBLIC WITH_QD)
 add_library(qd STATIC ${QdSources})
-target_include_directories(qd PUBLIC qd-2.3.22/include)
-target_include_directories(qd PUBLIC qd-2.3.22)
+target_include_directories(qd PUBLIC qd-2.3.22/include qd-2.3.22)
 
 #target_include_directories(mandel PUBLI#C qd-2.3.22/include)
 target_link_libraries(mandel PUBLIC qd)

+ 12 - 1
libmandel/include/ClGenerators.h

@@ -8,7 +8,7 @@
 #ifdef __APPLE__
 #include <OpenCL/cl.hpp>
 #else
-#include <CL/cl.hpp>
+#include <CL/cl2.hpp>
 #endif
 
 namespace mnd
@@ -61,6 +61,17 @@ protected:
 };
 
 
+class mnd::ClGeneratorDoubleDouble : public ClGenerator
+{
+public:
+    ClGeneratorDoubleDouble(cl::Device device, bool smooth);
+    virtual ~ClGeneratorDouble(void) = default;
+
+    virtual void generate(const MandelInfo& info, float* data);
+protected:
+    virtual std::string getKernelCode(bool smooth) const;
+};
+
 class mnd::ClGenerator128 : public ClGenerator
 {
 public:

+ 2 - 2
libmandel/include/Fixed.h

@@ -458,7 +458,7 @@ struct Fixed64
     ~Fixed64() = default;
 
 
-    inline Fixed64(uint64_t bits, bool dummy) :
+    inline Fixed64(uint64_t bits, bool /* dummy */) :
         bits{ bits }
     {
     }
@@ -502,7 +502,7 @@ struct Fixed64
         int32_t newLo = loup & 0xFFFFFFFF + (lolo >> 32);*/
         double d = int32_t(bits >> 32) + double(uint32_t(bits)) / (1ULL << 32);
         double od = int32_t(other.bits >> 32) + double(uint32_t(other.bits)) / (1ULL << 32);
-        return d * od * (other.sign != sign) ? -1 : 1;
+        return d * od * ((other.sign != sign) ? -1 : 1);
 
         //return Fixed64{ (uint64_t(newUp) << 32) | newLo, true };
     }

+ 1 - 0
libmandel/include/Generators.h

@@ -4,6 +4,7 @@
 #include "MandelUtil.h"
 
 #include <vector>
+#include <map>
 #include <utility>
 
 

+ 1 - 0
libmandel/include/Types.h

@@ -12,6 +12,7 @@
 //#       include <boost/multiprecision/float128.hpp>
 #   endif
 #   include <boost/multiprecision/cpp_int.hpp>
+#   include <boost/functional/hash.hpp>
 #endif
 
 #ifdef WITH_QD

+ 74 - 0
libmandel/src/ClGenerators.cpp

@@ -305,6 +305,80 @@ std::string ClGeneratorDouble::getKernelCode(bool smooth) const
 }
 
 
+ClGeneratorDoubleDouble::ClGeneratorDoubleDouble(cl::Device device, bool smooth) :
+    ClGenerator{ device }
+{
+    context = Context{ device };
+    Program::Sources sources;
+
+    std::string kcode = this->getKernelCode(smooth);
+
+    sources.push_back({ kcode.c_str(), kcode.length() });
+
+    program = Program{ context, sources };
+    if (program.build({ device }) != CL_SUCCESS) {
+        throw std::string(program.getBuildInfo<CL_PROGRAM_BUILD_LOG>(device));
+    }
+
+    queue = CommandQueue(context, device);
+}
+
+
+void ClGeneratorDouble::generate(const mnd::MandelInfo& info, float* data)
+{
+    ::size_t bufferSize = info.bWidth * info.bHeight * sizeof(float);
+
+    Buffer buffer_A(context, CL_MEM_WRITE_ONLY, bufferSize);
+    double pixelScaleX = double(info.view.width / info.bWidth);
+    double pixelScaleY = double(info.view.height / info.bHeight);
+
+    Kernel iterate = Kernel(program, "iterate");
+    iterate.setArg(0, buffer_A);
+    iterate.setArg(1, int(info.bWidth));
+    iterate.setArg(2, double(info.view.x));
+    iterate.setArg(3, double(info.view.y));
+    iterate.setArg(4, double(pixelScaleX));
+    iterate.setArg(5, double(pixelScaleY));
+    iterate.setArg(6, int(info.maxIter));
+
+    cl_int result = queue.enqueueNDRangeKernel(iterate, 0, NDRange(info.bWidth * info.bHeight));
+    queue.enqueueReadBuffer(buffer_A, CL_TRUE, 0, bufferSize, data);
+}
+
+
+std::string ClGeneratorDouble::getKernelCode(bool smooth) const
+{
+    return
+        "#pragma OPENCL EXTENSION cl_khr_fp64 : enable\n"
+        "__kernel void iterate(__global float* A, const int width, double xl, double yt, double pixelScaleX, double pixelScaleY, int max) {\n"
+        "   int index = get_global_id(0);\n"
+        "   int x = index % width;"
+        "   int y = index / width;"
+        "   double a = x * pixelScaleX + xl;"
+        "   double b = y * pixelScaleY + yt;"
+        "   double ca = a;"
+        "   double cb = b;"
+        ""
+        "   int n = 0;"
+        "   while (n < max - 1) {"
+        "       double aa = a * a;"
+        "       double bb = b * b;"
+        "       double ab = a * b;"
+        "       if (aa + bb > 16) break;"
+        "       a = aa - bb + ca;"
+        "       b = 2 * ab + cb;"
+        "       n++;"
+        "   }\n"
+        // N + 1 - log (log  |Z(N)|) / log 2
+        "   if (n >= max - 1)\n"
+        "       A[index] = max;\n"
+        "   else"
+        "       A[index] = ((float)n);\n"
+        //            "   A[index] = ((float)n) + 1 - (a * a + b * b - 16) / (256 - 16);\n"
+        //        "   A[get_global_id(0)] = 5;"
+        "}";
+}
+
 
 ClGenerator128::ClGenerator128(cl::Device device, bool smooth) :
     ClGenerator{ device }

+ 10 - 7
libmandel/src/Mandel.cpp

@@ -126,11 +126,6 @@ MandelContext::MandelContext(void) :
     adaptiveGenerator = std::make_unique<AdaptiveGenerator>();
     adaptiveGeneratorSmooth = std::make_unique<AdaptiveGenerator>();
 
-    adaptiveGenerator->addGenerator(1.0e-7, *cpuGeneratorFloat);
-    adaptiveGenerator->addGenerator(0.5e-15, *cpuGeneratorDouble);
-    adaptiveGeneratorSmooth->addGenerator(1.0e-7, *cpuGeneratorFloatSmooth);
-    adaptiveGeneratorSmooth->addGenerator(0.5e-15, *cpuGeneratorDoubleSmooth);
-
     {
         auto& device1 = devices[0];
         Generator* floatGenerator = device1.getGeneratorFloat(false);
@@ -139,19 +134,27 @@ MandelContext::MandelContext(void) :
         Generator* doubleGeneratorSmooth = device1.getGeneratorDouble(true);
         if (floatGenerator != nullptr)
             adaptiveGenerator->addGenerator(1.0e-7, *floatGenerator);
+        else
+            adaptiveGenerator->addGenerator(1.0e-7, *cpuGeneratorFloat);
         if (doubleGenerator != nullptr)
             adaptiveGenerator->addGenerator(0.5e-15, *doubleGenerator);
+        else
+            adaptiveGenerator->addGenerator(0.5e-15, *cpuGeneratorDouble);
         if (floatGeneratorSmooth != nullptr)
             adaptiveGeneratorSmooth->addGenerator(1.0e-7, *floatGeneratorSmooth);
+        else
+            adaptiveGeneratorSmooth->addGenerator(1.0e-7, *cpuGeneratorFloatSmooth);
         if (doubleGeneratorSmooth != nullptr)
             adaptiveGeneratorSmooth->addGenerator(0.5e-15, *doubleGeneratorSmooth);
+        else
+            adaptiveGeneratorSmooth->addGenerator(0.5e-15, *cpuGeneratorDoubleSmooth);
     }
 
 #ifdef WITH_QD
         adaptiveGenerator->addGenerator(Real("1.0e-29"), *cpuGeneratorDD);
         adaptiveGeneratorSmooth->addGenerator(Real("1.0e-29"), *cpuGeneratorDDSmooth);
-        adaptiveGenerator->addGenerator(Real("1.0e-50"), *cpuGeneratorQD);
-        adaptiveGeneratorSmooth->addGenerator(Real("1.0e-50"), *cpuGeneratorQDSmooth);
+        adaptiveGenerator->addGenerator(Real("1.0e-57"), *cpuGeneratorQD);
+        adaptiveGeneratorSmooth->addGenerator(Real("1.0e-57"), *cpuGeneratorQDSmooth);
 #endif
 #ifdef WITH_BOOST
         //adaptiveGenerator->addGenerator(1.0e-28, *cpuGeneratorQuad);

+ 85 - 0
libmandel/src/doubledouble.cl

@@ -0,0 +1,85 @@
+#pragma OPENCL EXTENSION cl_khr_fp64 : enable
+
+inline double2 twoSum(double a, double b) {
+  double s = a + b;
+  double bb = s - a;
+  double e = (a - (s - bb)) + (b - bb);
+  return (s, e);
+}
+
+inline double2 quickTwoSum(double a, double b) {
+  double s = a + b;
+  double e = b - (s - a);
+  return (s, e);
+}
+
+inline double2 twoProd(double a, double b) {
+//#ifdef QD_FMS
+  double p = a * b;
+  double e = fma(a, b, -p);
+  return (p, e);
+//#else
+//  double a_hi, a_lo, b_hi, b_lo;
+//  double p = a * b;
+//  split(a, a_hi, a_lo);
+//  split(b, b_hi, b_lo);
+//  err = ((a_hi * b_hi - p) + a_hi * b_lo + a_lo * b_hi) + a_lo * b_lo;
+//  return p;
+//#endif
+}
+
+inline double2 mul(double2 a, double2 b) {
+    double2 p = twoProd(a.0, b.0);
+    p.1 += a.0 * b.1 + a.1 * b.0;
+    return quickTwoSum(p.0, p.1);
+}
+
+inline double2 add(double2 a, double2 b) {
+    double se = twoSum(a.0, b.0);
+    se.1 += a.1 + b.1;
+    return quickTwoSum(se.0, se.1);
+}
+
+inline double2 mulDouble(double2 a, double b) {
+    double2 p = twoProd(a.0, b);
+    p.1 += a.1 * b;
+    return quickTwoSum(p.0, p.1);
+}
+
+__kernel void iterate(__global float* A, const int width,
+                      double x1, double x2, double y1, double y2,
+                      double pw1, double pw2, double ph1, double ph2, int max) {
+   int index = get_global_id(0);
+   int px = index % width;
+   int py = index / width;
+
+   double2 xl = (x1, x2);
+   double2 yt = (y1, y2);
+   double2 pixelScaleX = (pw1, pw2);
+   double2 pixelScaleY = (ph1, ph2);
+
+   double2 a = add(mulDouble(pixelScaleX, (double) px), xl); // pixelScaleX * px + xl
+   double2 b = add(mulDouble(pixelScaleY, (double) py), yt) // pixelScaleY * py + yt
+   double2 ca = a;
+   double2 cb = b;
+
+   int n = 0;
+   while (n < max - 1) {
+       double2 aa = mul(a, a);
+       double2 bb = mul(b, b);
+       double2 ab = mul(a, b);
+       if (aa.0 + aa.1 + bb.0 + bb.1 > 16) break;
+       double2 minusbb = (-bb.0, -bb.1)
+       a = add(add(aa, minusbb), ca);
+       double2 halfb = add(ab + cb);
+       b = add(halfb, halfb);
+       n++;
+   }
+// N + 1 - log (log  |Z(N)|) / log 2
+   if (n >= max - 1)
+       A[index] = max;
+   else
+       A[index] = ((float)n);
+//               A[index] = ((float)n) + 1 - (a * a + b * b - 16) / (256 - 16);
+//           A[get_global_id(0)] = 5;
+};