这篇文章是吧很久以前学的一个神奇算法归一下档，在公交车上突然想起来了，觉得还是很有必要再仔细梳理一遍，对以后也许有用。
先看图再说话：
Gray Scott Reaction-Diffusion算法， 在模拟微观细胞的运动或者类似的效果是非常神奇。
理论链接：http://www.karlsims.com/rd.html
原理：模拟两种物质之间在平面（暂时是平面）上的相互作用，动作分为反应与扩散。
公式：
右手边第一项为扩散，第二项为反应，第三项是供给，DA与DB为分散率，▽2A或B为相邻位置扩散阵列（下图），f和k是两个非常重要的外参数，分别叫做feed和kill。可以简单理解成A是环境中的实物，两个B和一个A就能变成三个B，过程中feed会变出更多的A，k也会去掉一些B。非常想现实的生长环境。
|0.05| 0.2 |0.05|
| 0.2 | -1 | 0.2 |
|0.05| 0.2 |0.05|
上面是拉普拉斯方程用到的阵列，因为外围总和等于中间数的反数，所以参数化后便可以随便调
在houdini里面，使用sop solver来迭代A和B的物质变化，在进入solver前，先定义好A和B的值，一般是A在整个平面是1.0，B随机散布开，这样出来的效果也好,如下图是我做的B随机散布示意图：
同时用noise方法随机分布feed和kill的值，你当然可以使用固定的f和k，不过有随机变化的话会让生长更加自然有趣。下面是detail面板里面的简单示意：
值得注意的是feed和kill的值是非常敏感的，我自己测试时，也就几个非常窄的区间内会有生长效果，要不然cb会非常快的全死掉。
在solver里面给preframe 加一个新的vex node，代码如下：
#pragma label resx "Res X"
#pragma label resy "Res Y"
#pragma label diff_a "Diffusion A"
#pragma label diff_b "Diffusion B"
#pragma label weight "Weight"
#pragma range weight 0.125 0.25
#pragma label feed "Feed"
#pragma range feed 0 1
#pragma label kill "Kill"
#pragma range kill 0 1
#pragma label time_step "Sub Step"
#pragma range time_step 0 20
#pragma hint ca invisible
#pragma hint cb invisible
//this fuction will diffuse
float laplacian(int ptnum; int resx; int resy; string attrib; float weight){
// | w_s | w | w_s |
//-----------------------
// | w |-1 | w |
//-----------------------
// | w_s | w | w_s | weight
// | pt+w-1 | pt+w | pt+w+1 |
//-----------------------------
// | pt-1 | pt | pt+1 |
//-----------------------------
// | pt-w-1 | pt-w | pt-w+1 | location
int vb = resx*2, hb = 1;
// detects boundaries
if((ptnum % resx) == 0 || (ptnum % resx) == (resx - 1)){
hb *= -1; // horizon
}
if(ptnum < resx || ptnum > (resy-1)*resx){
vb *= -1; // vetical
}
float n[];
resize(n,9);
float weight_second = 0.25 - weight;
n[0] = point(0, attrib, ptnum - vb - 1) * weight_second; // left bottom
n[1] = point(0, attrib, ptnum - vb ) * weight; // bottom
n[2] = point(0, attrib, ptnum - vb + 1) * weight_second; // right bottom
n[3] = point(0, attrib, ptnum - 1) * weight; // left
n[4] = point(0, attrib, ptnum ) * (-1.0); // center
n[5] = point(0, attrib, ptnum + 1) * weight; // right
n[6] = point(0, attrib, ptnum + vb - 1) * weight_second; // left top
n[7] = point(0, attrib, ptnum + vb ) * weight; // top
n[8] = point(0, attrib, ptnum + vb + 1) * weight_second; // right top
float sum=0.0;
foreach(float i; n){
sum += i;
}
return sum;
}
sop
gray_scott(
int resx = 512;
int resy = 512;
float diff_a = 1.0;
float diff_b = 1.0;
float weight = 0.15;
export float feed = 0.05;
export float kill = 0.05;
float time_step = 1.0;
export float ca = 0.0;
export float cb = 0.0;
)
{
float react_rate = ca * pow(cb,2);
float feed_set = feed * (1.0 - ca);
float kill_set = (kill + feed) * cb;
//float substep = 1.0/(float)time_step;
ca += (diff_a * laplacian(ptnum, resx, resy, "ca", weight) - react_rate + feed_set) * time_step;
cb += (diff_b * laplacian(ptnum, resx, resy, "cb", weight) + react_rate - kill_set) * time_step;
}
最后看一看不同的feed和kill值产生的不同效果：
feed = 0.051 kill = 0.062
feed = 0.0367 kill = 0.069
feed = 0.0561 kill = 0.0636
Houdini 中 Gray Scott Reaction-Diffusion算法的实现
原文：http://www.cnblogs.com/simonxia/p/4282569.html