gpu学习算法：生长的力量GPU
gpu学习算法：生长的力量GPU（来自百度百科）简单的说，分形就是研究无限复杂具备自相似结构的几何学。是大自然复杂表面下的内在数学秩序。我的 GPU 是 NVIDIA®Titan V，安置在一个比笔记本电脑更大的独立外壳中，并提供独立电源和散热装置。计算这些分形的速度大约是 CPU 的 300 倍，完全改变了我与 mandelbrot 程序的交互方式。我在程序中引入了 gpuArray，现在可以使用与屏幕分辨率相当的网格分辨率，并不断迭代，直到图像中的细节变得清晰可见。分形几何学（分形学）分形几何学是一门以不规则几何形态为研究对象的几何学。相对于传统几何学的研究对象为整数维数，如，零维的点、一维的线、二维的面、三维的立体乃至四维的时空。分形几何学的研究对象为非负实数维数，如0.63、1.58、2.72、log2/log3（参见康托尔集）。因为它的研究对象普遍存在于自然界中，因此分形几何学又被称为“大自然的几何学”。
Cleve Moler
MathWorks 的联合创始人兼首席科学家。Moler 曾是一名数学和计算机科学领域的教授，在密歇根大学、斯坦福大学和新墨西哥大学执教将近 20 年。除了是 MATLAB 第一版的作者之外，他还是 LINPACK 和 EISPACK 科学子程序库的作者之一。他和别人合著了三本有关数值方法的教科书。
我是个分形迷。
现在，我有了一个让我上瘾的 GPU。GPU 的设计初衷是图形加速，但 MATLAB 用它们来加速计算。分形学是需要大量计算的图形学，可以与 GPU 实现完美结合。
我的 GPU 是 NVIDIA®Titan V，安置在一个比笔记本电脑更大的独立外壳中，并提供独立电源和散热装置。计算这些分形的速度大约是 CPU 的 300 倍，完全改变了我与 mandelbrot 程序的交互方式。我在程序中引入了 gpuArray，现在可以使用与屏幕分辨率相当的网格分辨率，并不断迭代，直到图像中的细节变得清晰可见。
分形几何学（分形学）
分形几何学是一门以不规则几何形态为研究对象的几何学。相对于传统几何学的研究对象为整数维数，如，零维的点、一维的线、二维的面、三维的立体乃至四维的时空。分形几何学的研究对象为非负实数维数，如0.63、1.58、2.72、log2/log3（参见康托尔集）。因为它的研究对象普遍存在于自然界中，因此分形几何学又被称为“大自然的几何学”。
简单的说，分形就是研究无限复杂具备自相似结构的几何学。是大自然复杂表面下的内在数学秩序。
（来自百度百科）
分形这个术语是由 Benoit Mandelbrot (1924-2010) 提出的，他是一位波兰/法国/美国的数学家，职业生涯的大部分时间都在位于纽约约克敦海茨的 IBM 华生研究中心度过。
GPU 重新激起了我对曼德尔勃特集合本身的兴趣，并通过三种分形变体简化了我的工作，这三种变体包括：“燃烧之船”、“权力之塔”和牛顿法的全局收敛行为。
◆ ◆ ◆ ◆
曼德尔勃特集合
1978 年，Robert W. Brooks 和 Peter Matelsk 在一台行式打印机上打印出了曼德尔勃特集合的首张图片。
曼德尔勃特集合的行式打印机图像。
第一张彩色图片出现在 1985 年的《科学美国人》封面上，那时计算机图形显示正在普及。
曼德尔勃特集合的彩色图像。
从那时起，曼德尔勃特集合激发了数学爱好者对数学的深入研究，并成为许多图形项目、硬件演示和 Web 页面的基础案例。
曼德尔勃特集合解密
曼德尔勃特现象一直存在。
曼德尔勃特涉及到复数的简单迭代，从初始点 z0 开始。曼德尔勃特集合是由 z0 值组成的复平面上的区域，其中由以下公式定义的轨迹保持有界：
这就是整个定义。如此简单的定义竟能产生如此迷人的复杂性，真是令人惊讶。
上面所示的彩色图像是一个曼德尔勃特集合的鸟瞰视图。它没有办法显示集合边界外缘详细的结构。事实上，曼德尔勃特集合有小的细丝达到边缘地区，但在此视图中很难看到边缘。
图 1 是曼德尔勃特集合几何图形的略图。
最大的组成部分是一个心形的心脏曲线，由以下参数方程的曲线包围：
图1. 曼德尔勃特集合几何图形的略图。
心脏线的左边是一个半径为 ¼ 的圆盘。心脏线和圆盘的连接部位在红点处。稍后再详细介绍那个红点。
在这两个部件的周围是无穷多的近似圆形的圆盘，圆盘的直径逐渐递减。最小圆盘的精确位置和形状只能通过计算来确定。
最近已经证明，曼德尔勃特集合是数学上连接的，但连接区域有时太薄，我们无法在图形屏幕上解析它，甚至无法在合理的时间内计算它。
计算曼德尔勃特集合
我们联想到的曼德尔勃特的迷人图案来自集合外的边缘区域。
这是计算的核心函数。输入是一个复标量起点 z0 和一个实标量参数，我称之为 depth。输出是一个计数 kz。如果因为 z 在半径为 2 的圆盘之外而终止迭代，则 z 注定会变成无穷大，而计数 kz 可以用作色图的索引。另一方面，如果 z 在 depth 迭代中幸存，则声明它在曼德尔勃特集合中。
function kz = mandelbrot(z0 depth) z = z0; kz = 0; while (z*conj(z) <= 4) & (kz <= depth) kz = kz 1; z = z*z z0; end end
让我们在复平面上生成一个点网格，覆盖宽度为 w 的平方，以复点 zc 为中心。
grid = 512; s = w*(-1:1/grid:1); [u v] = meshgrid(s real(zc) s imag(zc));
现在要介绍的是，在 CPU 上运行的程序和在 GPU 上运行的程序之间的唯一区别。要在 CPU 上生成网格：
z0 = u i*v;
对于GPU：
z0 = gpuArray(u i*v);
然后，对于任一处理器，以下语句
kz = arrayfun(@mandelbrot z0);
将标量函数 mandelbrot 应用于数组 z0 的所有元素。在 CPU 上，这是一个向量化网格上的双重 for 循环。在 GPU 上，语句被分解成数百个并行运行的独立任务。
现在，让我们来看看来自曼德尔勃特边缘的两幅不同寻常的图像。
在边缘上
被称为海马谷的区域位于中央心脏线和它左边的圆盘之间。实际上有两个山谷，一个在实轴上，一个在实轴下。它们在图 1 中看到的红点处交汇。发挥点想象力，您就可以想象出生活在图 2 中的小海洋生物。稍后我们将看到，山谷也包含隐伏的 π。
图2. 海马谷。
由于曼德尔勃特集合是自相似的，所以它包含无限多个微型曼德尔勃特，每个曼德尔勃特的形状都与大的曼德尔勃特相同。如图 3 所示，其放大倍数为 1010。
图3. 微型曼德尔勃特。
隐伏的 π
1991 年，科罗拉多州立大学的 Dave Boll 发现了一个惊人的结果。Boll 当时正在研究曼德尔勃特迭代在海马谷中的行为。山谷接近轴时会变得更窄，它们在图 1 中的红点处 (-3/4 0) 交汇。
我们可以在一个以 -3/4 yi 为中心、刚刚偏离轴的小网格上重复 Boll 的计算，y 代表一个极小的值，i代表虚数单位：
y = 1.0e-07 zc = -3/4 y*i
使网格足够大，刚好可以接触到轴：
width = 2*y grid = 4
选择 depth 使其与 y 成反比，这会使它变得非常大：
depth = 4/y
使用这些参数，运行上面的代码。它会生成：
kz = 40000000 40000000 40000000 40000000 40000000 40000000 40000000 40000000 40000000 40000000 40000000 40000000 31415926 40000000 40000000 40000000 40000000 20943951 40000000 40000000 40000000 40000000 15707963 40000000 40000000
查看网格中心的迭代计数。是否看到了一个熟悉的值？
这不是侥幸。Aaron Klebanof 于 2001 年发表的一篇论文《曼德尔勃特集合中的 π》分析了心脏线前端尖点的类似计算。
接下来，我们要介绍的是曼德尔勃特迭代的一个奇怪变体。
燃烧之船
燃烧之船来自一个奇怪的迭代：
其中
我说这很奇怪，是因为函数 F(z) 不是解析的。我对这个迭代很感兴趣，因为在以下图片中有不可思议的相似之处。
如图 4 所示，初始域是一个以 -0.5-0.5i 为中心的宽度为 3.5 的正方形。我在船的尾迹处插入了一个指向相关区域的箭头。
图4. 燃烧的船，初始域。
将船的尾迹放大 500 倍至 -1.861-.002i。使用 bone 色图使它看起来是冷的而不是燃烧的。图 5 显示了分形的结果，旁边是 1915 年南极探险家欧内斯特·沙克尔顿的船“耐力”号在威德尔海冻结的照片。
图5. 左图：放大后显示的燃烧之船的尾迹。右图：Hurley 于 1915 年拍摄的沙克尔顿的船在南极洲结冰的照片。
图片来源：澳大利亚国家图书馆 http://nla.gov.au/nla.obj-158931586/view
权利之塔
从任意复数 z 开始，反复取幂。
我们可以将其表示为迭代。从 y0=1 开始，令
如果 z 太大，则此迭代将膨胀到无穷大。对于某些 z，它会收敛到一个有限的极值。例如，如果z=√2，yk 将收敛于 2。最有趣的情况是 yk 接近一个循环。例如，如果 z 接近 2.5 4i，则循环长度为 7。
2.4684 4.0754i -0.6216 0.3634i 0.2603 - 0.0184i 1.4868 0.3613i -3.4877 6.1054i 7.7632e-06 - 2.6617e-06i 1.0000 0.0000i 2.4684 4.0755i
这个循环长度是“权力之塔”分形的基础。图 6 显示了整体的分形。
图6. 权力之塔分形。
放大 105 倍并更改色图将生成图 7 所示的图像。
图7. 权力之塔分形的细节。
牛顿法
当牛顿法的起始点不接近函数的零点时，全局行为会显得不可预测。从复杂平面上一个起始点区域到收敛点的迭代计数的等高线图生成了令人深思的分形图像。
迭代是常见的。
选择一个函数 f(z) 对 f′(z) 求导。从 z0 开始，令：
它最终会收敛到 f 的零点。计算接近目标所需的迭代次数。
有许多函数（和色图）可供选择。我最喜欢的三次多项式是：
图 8 显示了被划分为三个区域的复杂平面，其中迭代收敛于立方的三个零中的一个。这些区域之间是强烈的分形作用区域，如图中黑色部分所示。
图8. 在 z3– 2z – 5 上的牛顿的迭代。
图 9 显示了寻找以下函数式零点的牛顿法的全局行为：
该函数有无穷多个零和一个无界的一阶导数。
图9. 有关tan(sin(z)) – sin(tan(z))牛顿迭代的详细信息。
图 10 的函数是：
最突出的蓝色区域环绕在 ±1/π 的零处。
图10. 关于 z sin(1/z) 的牛顿迭代的详细信息。
本文所示的许多图像对我来说是全新的——我以前从未见过它们。利用 GPU 进行计算，使它们得以展示在公众面前。
◆ ◆ ◆ ◆
