当前位置: 首页 > 知识库问答 >
问题:

绘制mandelbrot集的最快算法是什么?

琴镜
2023-03-14

我已经尝试了很多算法来渲染Mandelbrot集,包括简单的逃逸时间算法,以及优化的逃逸时间算法。但是,有没有更快的算法可以像我们在YouTube上看到的那样有效地产生真正深的缩放。此外,我很想得到一些想法,如何提高我的精度超过C/C

共有3个答案

蒋奇
2023-03-14

我的最快解决方案通过遵循轮廓边界和填充避免在相同深度的大区域上迭代。有一个惩罚,它是有可能掐断小芽,而不是去他们周围,但总的来说,一个小的代价来支付快速缩放。

一种可能的效率是,如果缩放比例加倍,则已经有¼个点。

对于动画,我将每个帧的值归档,每次倍增比例,并在实时播放时插值帧间的值,因此动画每秒倍增一次。double类型允许存储50多个关键帧,使动画持续一分钟以上(先入后出)。

实际的迭代由手工制作的汇编程序完成,因此一个像素完全在FPU中迭代。

农均
2023-03-14

优化的转义算法应该足够快,以实时绘制Mandelbrot集。您可以使用多个线程,以便更快地实现(例如,使用OpenMP非常容易)。您还可以使用SIMD指令手动对代码进行矢量化,以便在需要时更快。您甚至可以使用着色器和/或GPU计算框架(OpenCL或CUDA)直接在GPU上运行此操作,如果这对您来说还不够快的话(尽管要有效地执行此操作有点复杂)。最后,您应该调整迭代次数,使其非常小。

缩放不应对性能产生任何直接影响。它只是改变了计算的输入窗口。但是,它确实会产生间接影响,因为实际迭代次数会发生变化。不应计算窗口外的点。

双精度也应足以正确绘制Mandelbrot集。但是如果你真的想要更精确的计算,你可以使用双精度,它可以提供相当好的精度,性能也不会太差。但是,手动实现双精度有点棘手,而且仍然比仅使用双精度慢得多。

端木鹏
2023-03-14

与普通GPU相比,即使是高端CPU也会慢得多。即使在GPU上使用简单的迭代算法,您也可以实时渲染。因此,在GPU上使用更好的算法可以获得高缩放,但是对于您需要的任何体面算法:

  • 多路渲染,因为我们不能在GPU上自行修改纹理
  • 高精度浮点数,因为浮点数/双倍数是不够的。

这里有几个相关的QA:

  • GLSL RT Mandelbrot

这可能会让你开始。。。

一种加速的方法是你可以像我在第一个链接中那样使用分数逃逸。它提高了图像质量,同时保持最大迭代次数较低。

第二个链接将使您大致了解分形的哪些部分在内部和外部,以及距离有多远。它不是很精确,但可以用来避免计算“肯定在外面”的零件的迭代。

下一个链接将向您展示如何实现更好的精度。

最后一个链接是关于微扰的,其思想是,您仅对一些参考点使用高精度数学,并使用低精度数学计算其相邻点,而不会失去精度。但从未使用过,但看起来很有希望。

最后,一旦你实现了快速渲染,你可能会想要实现这个目标:

  • 如何在缩放Mandelbrot集时调整平移

这里是GLSL中用于单值的3*64位双精度的一个小示例:

// high precision float (very slow)
dvec3 fnor(dvec3 a)
    {
    dvec3 c=a;
    if (abs(c.x)>1e-5){ c.y+=c.x; c.x=0.0; }
    if (abs(c.y)>1e+5){ c.z+=c.y; c.y=0.0; }
    return c;
    }
double fget(dvec3 a){ return a.x+a.y+a.z; }
dvec3 fset(double a){ return fnor(dvec3(a,0.0,0.0)); }
dvec3 fadd(dvec3 a,double b){ return fnor(a+fset(b)); }
dvec3 fsub(dvec3 a,double b){ return fnor(a-fset(b)); }
dvec3 fmul(dvec3 a,double b){ return fnor(a*b); }
dvec3 fadd(dvec3 a,dvec3 b){ return fnor(a+b); }
dvec3 fsub(dvec3 a,dvec3 b){ return fnor(a-b); }
dvec3 fmul(dvec3 a,dvec3 b)
    {
    dvec3 c;
    c =fnor(a*b.x);
    c+=fnor(a*b.y);
    c+=fnor(a*b.z);
    return fnor(c);
    }

所以每个hi精度值是dvec3... fnor中的阈值可以更改为任何范围。您可以将其转换为vec3浮动...

[Edit1]“快速”C示例

好的,我想尝试我的新SSD1306驱动程序和我的AVR32 MCU来计算Mandelbrot,这样我就可以比较速度和这个Arduino 3D Pong Mandelbrot。我使用AT32UC3A3256,带有~66MHz无FPU无GPU和128x64x1bpp显示器。没有外部存储器,只有内部16 32 KB。NaiveMandlebrot的速度很慢(每帧约2.5秒),所以我做了这样的事情(利用这个位置,视图的缩放是连续的):

>

抖动留出空间,因为我的输出仅为B

使用变量max迭代n基于zoom

在更改n时,使最后一帧无效以强制执行完全重新计算。我知道这很慢,但它只发生3次变焦范围之间的转换。

从最后一帧的缩放计数看起来不太好,因为它不是线性的。

可以使用最后的计数,但也需要记住用于迭代的复杂变量,这会占用太多内存。

记住最后一帧,也记住哪个x, y屏幕坐标映射到哪个Mandelbrot坐标。

在每个帧上计算屏幕坐标和Mandelbrot坐标之间的映射。

重新映射最后一帧以调整到新的位置和缩放

因此,只需查看#3和#4中的数据,如果我们在最后一帧和实际帧中的位置相同(接近像素大小的一半),则复制像素。然后重新计算其余的。

如果您的视图平滑,这将极大地提高性能(因此位置和缩放不会在每帧基础上发生很大变化)。

我知道它的描述有点模糊,所以这里有一个C代码,你可以推断出所有的疑问:

//---------------------------------------------------------------------------
//--- Fast Mandelbrot set ver: 1.000 ----------------------------------------
//---------------------------------------------------------------------------
template<int xs,int ys,int sh> void mandelbrot_draw(float mx,float my,float zoom)
    {
    // xs,ys - screen resolution
    // sh    - log2(pixel_size) ... dithering pixel size
    // mx,my - Mandelbrot position (center of view) <-1.5,+0.5>,<-1.0,+1.0>
    // zoom  - zoom
    // ----------------
    // (previous/actual) frame
    static U8 p[xs>>sh][ys>>sh];            // intensities (raw Mandelbrot image)
    static int n0=0;                        // max iteraions
    static float px[(xs>>sh)+1]={-1000.0};  // pixel x position in Mandlebrot
    static float py[(ys>>sh)+1];            // pixel y position in Mandlebrot
    // temp variables
    U8 shd;                                 // just pattern for dithering
    int ix,iy,i,n,jx,jy,kx,ky,sz;           // index variables
    int nx=xs>>sh,ny=ys>>sh;                // real Mandelbrot resolution
    float fx,fy,fd;                         // floating Mandlebrot position and pixel step
    float x,y,xx,yy,q;                      // Mandelbrot iteration stuff (this need to be high precision)
    int   qx[xs>>sh],qy[ys>>sh];            // maping of pixels between last and actual frame
    float px0[xs>>sh],py0[ys>>sh];          // pixel position in Mandlebrot from last frame
    // init vars
         if (zoom<   10.0) n= 31;
    else if (zoom<  100.0) n= 63;
    else if (zoom< 1000.0) n=127;
    else                   n=255;
    sz=1<<sh;
    ix=xs; if (ix>ys) ix=ys; ix/=sz;
    fd=2.0/(float(ix-1)*zoom);
    mx-=float(xs>>(1+sh))*fd;
    my-=float(ys>>(1+sh))*fd;
    // init buffers
    if ((px[0]<-999.0)||(n0!=n))
        {
        n0=n;
        for (ix=0;ix<nx;ix++) px[ix]=-999.0; 
        for (iy=0;iy<ny;iy++) py[iy]=-999.0;
        for (ix=0;ix<nx;ix++)
         for (iy=0;iy<ny;iy++)
          p[ix][iy]=0;
        }
    // store old and compute new float positions of pixels in Mandelbrot to px[],py[],px0[],py0[]
    for (fx=mx,ix=0;ix<nx;ix++,fx+=fd){ px0[ix]=px[ix]; px[ix]=fx; qx[ix]=-1; }
    for (fy=my,iy=0;iy<ny;iy++,fy+=fd){ py0[iy]=py[iy]; py[iy]=fy; qy[iy]=-1; }
    // match old and new x coordinates to qx[]
    for (ix=0,jx=0;(ix<nx)&&(jx<nx);)
        {
        x=px[ix]; y=px0[jx];
        xx=(x-y)/fd; if (xx<0.0) xx=-xx;
        if (xx<=0.5){ qx[ix]=jx; px[ix]=y; }
        if (x<y) ix++; else jx++;
        }
    // match old and new y coordinates to qy[]
    for (ix=0,jx=0;(ix<ny)&&(jx<ny);)
        {
        x=py[ix]; y=py0[jx];
        xx=(x-y)/fd; if (xx<0.0) xx=-xx;
        if (xx<=0.5){ qy[ix]=jx; py[ix]=y; }
        if (x<y) ix++; else jx++;
        }
    // remap p[][] by qx[]
    for (ix=0,jx=nx-1;ix<nx;ix++,jx--)
        {
        i=qx[ix]; if ((i>=0)&&(i>=ix)) for (iy=0;iy<ny;iy++) p[ix][iy]=p[i][iy];
        i=qx[jx]; if ((i>=0)&&(i<=jx)) for (iy=0;iy<ny;iy++) p[jx][iy]=p[i][iy];
        }
    // remap p[][] by qy[]
    for (iy=0,jy=ny-1;iy<ny;iy++,jy--)
        {
        i=qy[iy]; if ((i>=0)&&(i>=iy)) for (ix=0;ix<nx;ix++) p[ix][iy]=p[ix][i];
        i=qy[jy]; if ((i>=0)&&(i<=jy)) for (ix=0;ix<nx;ix++) p[ix][jy]=p[ix][i];
        }
    // Mandelbrot
    for (iy=0,ky=0,fy=py[iy];iy<ny;iy++,ky+=sz,fy=py[iy]) if ((fy>=-1.0)&&(fy<=+1.0))
     for (ix=0,kx=0,fx=px[ix];ix<nx;ix++,kx+=sz,fx=px[ix]) if ((fx>=-1.5)&&(fx<=+0.5))
        {
        // invalid qx,qy ... recompute Mandelbrot
        if ((qx[ix]<0)||(qy[iy]<0))
            {
            for (x=0.0,y=0.0,xx=0.0,yy=0.0,i=0;(i<n)&&(xx+yy<4.0);i++)
                {
                q=xx-yy+fx;
                y=(2.0*x*y)+fy;
                x=q;
                xx=x*x;
                yy=y*y;
                }
            i=(16*i)/(n-1); if (i>16) i=16; if (i<0) i=0;
            i=16-i; p[ix][iy]=i;
            }
        // use stored intensity
        else i=p[ix][iy];
        // render point with intensity i coresponding to ix,iy position in map
        for (i<<=3                 ,jy=0;jy<sz;jy++)
         for (shd=shade8x8[i+(jy&7)],jx=0;jx<sz;jx++)
          lcd.pixel(kx+jx,ky+jy,shd&(1<<(jx&7)));
        }
    }
//---------------------------------------------------------------------------
//---------------------------------------------------------------------------
//---------------------------------------------------------------------------

可在链接的SSD1306 QA中找到lcdshade8x8内容。但是,您可以忽略它,它只是抖动并输出一个像素,因此您可以直接输出i(即使没有缩放到

此处预览(在PC上,因为我懒得连接摄像头…):

因此,它将64x32个Mandelbrot像素显示为128x64抖动图像。在我的AVR32上,这可能比简单方法快8倍(可能3-4fps)。。。代码可能会更优化,但请记住,Mandelbrot并不是唯一运行的东西,因为我有一些ISR处理程序在后台处理LCD,还有基于此的TTS引擎,从那时起我对其进行了大量升级,并用于调试(是的,它可以与渲染并行)。此外,我的3D引擎占用了大量~11 KB的内存(大部分是深度缓冲区),内存不足。

预览是使用以下代码完成的(内部计时器):

static float zoom=1.0;
mandelbrot_draw<128,64,1>(+0.37,-0.1,zoom);
zoom*=1.02; if (zoom>100000) zoom=1.0;

同样对于非AVR32 C环境使用这个:

//------------------------------------------------------------------------------------------
#ifndef _AVR32_compiler_h
#define _AVR32_compiler_h
//------------------------------------------------------------------------------------------
typedef int32_t  S32;
typedef int16_t  S16;
typedef int8_t   S8;
typedef uint32_t U32;
typedef uint16_t U16;
typedef uint8_t  U8;
//------------------------------------------------------------------------------------------
#endif
//------------------------------------------------------------------------------------------

[Edit2]GLSL中的浮点精度更高

Mandelbrot的主要问题是,它需要添加指数差非常大的数字。对于,-操作,我们需要对齐两个操作数的尾数,将它们添加为整数,并将其规格化回科学记数法。但是,如果指数差很大,则结果尾数需要的位数超过32位浮点所能容纳的位数,因此只保留24个最高有效位。这将创建导致像素化的舍入错误。如果查看二进制中的32位float,您将看到:

float a=1000.000,b=0.00001,c=a+b; 

//012345678901234567890123456789 ... just to easy count bits
a=1111101000b                                         // a=1000
b=         0.00000000000000001010011111000101101011b  // b=0.00000999999974737875
c=1111101000.00000000000000001010011111000101101011b  // not rounded result
c=1111101000.00000000000000b                          // c=1000 rounded to 24 bits of mantissa

现在的想法是扩大尾数位数。最简单的技巧是使用2个浮点数,而不是一个浮点数:

//012345678901234567890123 ... just to easy count bits
a=1111101000b                                         // a=1000
                       //012345678901234567890123 ... just to easy count     b=         0.0000000000000000101001111100010110101100b    // b=0.00000999999974737875
c=1111101000.00000000000000001010011111000101101011b  // not rounded result
c=1111101000.00000000000000b                          // c=1000 rounded to 24 
 +          .0000000000000000101001111100010110101100b
                           //012345678901234567890123 ... just to easy count bits

所以结果的一部分在一个浮动中,其余的在另一个浮动中...每个值的浮点数越多,尾数就越大。然而,在GLSL中,将大尾数精确划分为24个比特块是复杂而缓慢的(如果可能的话,由于GLSL的限制)。相反,我们可以为每个浮点数选择一些指数范围(就像上面的例子一样)。

所以在这个例子中,我们得到了3个浮点数(vec3)每一个(浮点数)值。每个坐标代表不同的范围:

       abs(x) <= 1e-5
1e-5 < abs(y) <= 1e+5
1e+5 < abs(z) 

value=(xyz)所以我们有3*24位尾数,但是范围并不完全匹配24位。为此,指数范围应除以:

log10(2^24)=7.2247198959355486851297334733878

而不是10。。。例如,类似这样的事情:

       abs(x) <= 1e-7
1e-7 < abs(y) <= 1e+0
1e+0 < abs(z) 

此外,必须选择范围,以便它们处理您使用的值范围,否则它将是免费的。所以如果你的数字是<代码>

注意一些舍入(但远低于本机浮点)仍然发生!!!

现在,对这些数字进行操作比对普通浮点数进行操作稍微复杂一点,因为您需要将每个值作为括号中所有组件的总和进行处理,以便:

(x0+y0+z0) + (x1+y1+z1) = (x0+x1 + y0+y1 + z0+z1)
(x0+y0+z0) - (x1+y1+z1) = (x0-x1 + y0-y1 + z0-z1)
(x0+y0+z0) * (x1+y1+z1) = x0*(x1+y1+z1) + y0*(x1+y1+z1) + z0*(x1+y1+z1)

不要忘记将值规格化回定义的范围。避免添加大小(abs)值,以避免x0 z0等。。。

[Edit3]新的win32演示CPU与GPU

  • win32 Mandelbrot演示64位浮点数

两个可执行文件都预设在同一位置,并在doubles开始四舍五入时缩放显示。当y轴在这个位置开始偏离10^9左右时,我不得不稍微升级计算px,py坐标的方法(对于其他位置,阈值可能仍然太大)

这里预览CPU与GPU的高变焦(n=1600):

CPU的RT GIF捕获(n=62,按比例缩小GIF 4x):

 类似资料:
  • 我正在渲染一个mandelbrot集,并且已经实现了一些平滑的着色,但是当仔细观察时,图片变得非常嘈杂。我想知道,为了达到更好的美感,改善我的色彩的最好方法是什么。使用直方图着色是否有助于去除粗糙的像素区域?这里是一个使用10000次迭代的分形渲染。 这是我现在生成和分配颜色的方式:

  • 我知道关于这件事,阿雷迪回答了很多问题。然而,我的略有不同。无论何时我们实现平滑着色算法,我都理解它。 其中n是逃逸迭代,2是z的幂,如果我没有弄错,z是该逃逸迭代处复数的模。然后,我们在颜色之间的线性插值中使用这个重整化的逃逸值来生成平滑的带状mandelbrot集。我已经看到了关于这个的其他问题的答案,我们通过HSB到RGB的转换来运行这个值,但是我仍然无法理解这将如何提供平滑的颜色渐变,以及

  • 我正在编写一个Java程序来显示我的入门编程类的Mandelbrot集。我相信我已经正确地设置了所有的数学,但是当我尝试绘制分形时,我得到的只是一种纯色。我已经测试了数学,它似乎应该是有效的。我搜索了一个多小时,但没有找到任何有用的东西。下面是我的复数类,并实际创建了Mandelbrot集:复数 曼德布罗特 我已经做了一些JUnit测试,上面的两个类似乎都可以工作。我的测试中可能有一个缺陷导致了疏

  • 问题内容: 我希望根据计算出的像素值绘制图像,以可视化某些数据。本质上,我希望采用彩色三元组的二维矩阵并将其渲染。 请注意,这不是图像处理,因为我既不变换现有图像也不做任何形式的全图像变换,它也不是矢量图形,因为我要渲染的图像没有预定的结构,我可能一次要产生一个像素的无定形斑点。 我现在需要渲染大约1kx1k像素的图像,但是可伸缩的东西会很有用。最终目标格式为PNG或任何其他无损格式。 目前,我一

  • 我试图画曼德布罗特集,其中的点是黑色的,其他的都是白色的。在这个初始版本中,我不希望能够放大,而只是创建一个静态图像。 我创建了一个ComplexNumber类,如下所示,用于处理平方运算和将复数相加。 这是我渲染GUI并实际计算Mandelbrot Set中的点的代码。 运行完这段代码后,我得到了下图。看起来曼德尔布罗特的布景有一点模糊,但随后被一吨黑色遮住了。我做错了什么? 更新的解决方案如下

  • 问题内容: 设和为两个集合。我正在寻找一种 非常 快速或优雅的方法来计算它们之间的设置差异(或,取决于您的偏好)。如标题所示,这两组存储和存储为Javascript数组。 笔记: 壁虎特技可以 我更喜欢本机函数(但是如果速度更快,我可以使用轻量级库) 我看过但未测试JS.Set(请参阅上一点) 编辑: 我注意到有关包含重复元素的集合的评论。当我说“设置”时,我指的是数学定义,这意味着(除其他外)它