Houdini 18 的Vellum目前使用的是更高级的PBD: XPBD(eXtended Position Based Dynamics) 。XPBD使用了2nd Order Integration ，传统的PBD使用的1nd Order Intergration 。想要查看的话，如下图（默认是2nd Order Integration)。
Vellum Solver中因为要对约束重复的运算（Constraint Iterations 设置迭代的次数），需要比较高的效率。 为了能够让并行运算，引入了Graph Corlor技术。Graph Color说白了， 就是定义一个属性，每个约束的属性和与它相邻的约束的属性值不同。
比如用一个Grid，用Vellum Contraints 得到 distant contrainsts ，然后用Graph Color节点进行计算 （结果默认写入在一个整型color属性上），最后用wrangle 显示颜色（相同的color属性颜色相同）， 如下图，能够看到，每一条边（primitive）都与相邻的边的颜色不同。 这样让颜色相同的contrainst在同一个Workset。 不同颜色，不同的Workset。这样在对constrait进行迭代计算时（是在opencl里面进行的），Opencl会根据Workset，同一个Workset才会一起平行计算。
当然，Vellum 中Graph Color的计算时在Vellum Solver内部，这里只是为了演示方便。
另个案例了解Vellum Solver 中 的 Graph Color
新建个一个Grid，选中两边的点， 给Pin Constraints，把两边定住。加 distant constraints ( Vellum Constraint 选 Distant along edges)，设置如下图。注意Stiffness设置为最大(10*1e+10)，这种材质现实中是没有弹性的， 但是我们又把Rest Length Scale设置为0.5，Vellum Solver会迫使它压缩（现实中这种材质的布料是不存在的）。然后把Vellum Solver的 Smooth Iteration 设置为0（这步很重要）。 结算。
结算后发现，有的边的长度缩短大致0.5倍，有的反而被拉长了。用Graph Color查看，发现拉伸边长的是同一个颜色，即作为同一个WorkSet在opencl中并行计算的。如下图：
这是因为黄色的边做一个同一个Workset首先在opencl中运算，红色和蓝色在黄色后面分别进行运算，这样导致错误都转嫁到了第一次计算的边上（这里的错误是指，设置了Rest Length Scale为0.5倍，红色和蓝色都正确压缩了，黄色反而拉长了，这是错误的， 不是我们想要的）。
Graph Color 的算法分析
下面分析下Graph Color节点里面的算法。新建一个Grid，Graph Color， Color（random from attributes）。
进入Graph Color内部， 主要分析最左边这路（Connectivity 为Primtive），Contraints运算的时候也是以Primtive为单位。
（1）create_color_attrib
创建 color (Integer)primitive属性
（2）find_npts
把当前primitive总共有多少邻居记录在属性 @__npts属性上。
i@__npts = len(primpoints(0, @primnum));
（3）color_high_valence （Run Over Detail ）
这步主要是把邻居多的primitive筛选出来（Graph Color 的参数Max Valence设置阈值，默认20），单独计算。因为进入opencl进行并行计算时，需要把邻居primitive写入在一个array数组里，如果相邻的邻居太多，会造成存储和计算时间的极大浪费。挑其中比较重要的代码部分进行分析：
下面的代码把筛选出来的符合条件的primitve（@__npts>maxvalence） 存入在prims数组里面。
int maxvalence = chi("maxvalence"); string grp = sprintf("\@__npts>=%d", maxvalence); int prims[] = expandprimgroup(0, grp);
然后对数组进行排序。npts数组记录prims数组中对应的prim的点的数量乘以-1 。argsort对npts数组排序，它返回的是一个indice的数组。
indice： 比如一个整型数组 int a[] = {0,4,2} 。a[1] = 4 ，中括号里面的1就是indice, 也就是数组下标。
对于一个数组int b[] = {-4，-10，-2，-8} 进行 argsort排序会得到一个数组下标的数组 {1,3,0,2}，最后用reorder函数得到最终的排序数组 {b[1],b[3],b[0],b[2]} 即{ -10，-8，-4，-2}。 （argsort 和reorder一般搭档起来用，来对数组排序）。
foreach(int i; int prim; prims) { npts[i] = -len(primpoints(0, prim)); } prims = reorder(prims, argsort(npts));
最终得到一个数组，根据primitive的点的数量由大到小排列。
核心的部分在下面这个foreeach部分
foreach(int i; int prim; prims) { int useNewColor = 1; // 声明新变量，是否使用新颜色，默认1，即使用新染色 int pts[] = primpoints(0, prim); for(int j = 0; j < ncolors; j++) // ncolors记录当前不同color的数量， 遍历所有color { int useColor = 1; // 声明新变量，能否用当前颜色，默认1， 可以使用 foreach(int pt; pts) // 遍历当前primitive的每个点 { if (mapping[j * npt + pt] != 0) //mapping是一个数组，npts是所有点的数量（不止当前primitive） { //比如npts= 14 ，当前primitive2个点，对于j=0,即color值为0 ，mapping[0*14+0]~mappint[0*14+(14-1)] useColor = 0; //对应的所有的点使用的color 0的情况，mapping值若为1即不等于0， 说明当前color值被邻居占用，不能使用 break; // 这是把useColor赋值0， 即当前color值不可用 } } if (useColor) //如果当前颜色可用 { colors[i] = j; // 把当前primitive 的color值记录到colors数组里 foreach(int pt; pts) // 对于当前primitive的所有点pt,color为j,对应mapping[j*npt+pt],设置为1，表示被当前primitive占用
{ mapping[j * npt + pt] = 1; // } useNewColor = 0; // 这种情况就不使用新颜色了 break; } } if (useNewColor) //说明没有可用的color,需要新增加一个color { colors[i] = ncolors; // 新的color值存到colors数组里，colors[i]是prims[i]对应的颜色值 ncolors++; // 总颜色数加1 resize(mapping, ncolors * npt); // mapping数组resize,空间比之前多了npt个 foreach(int pt; pts) { mapping[colors[i] * npt + pt] = 1; //对于新增的color[i],当前primitive的点pt, 把对应的mappint[colors[i]*npt+pt]赋值1， 表示被占用 } } }
这个节点完整的代码如下：
1 // Handle high valence prims in a serial pre-pass, 2 // else the find_neighbors wrangle below ends up creating 3 // extremely large neighbors arrays, which are expensive 4 // in memory and time. 5 int maxvalence = chi("maxvalence"); 6 string grp = sprintf("\@__npts>=%d", maxvalence); 7 int prims[] = expandprimgroup(0, grp); 8 int nprim = len(prims); 9 if (nprim == 0) 10 return; 11 12 // Create mapping for each point to each color and 13 // prims to each color. 14 int npt = npoints(0); 15 int mapping[]; 16 int colors[]; 17 resize(colors, nprim); 18 int ncolors = 0; 19 20 21 // Reverse sort by num pts in constraint 22 // so constraints with more points are first. 23 int npts[]; 24 resize(npts, nprim); 25 foreach(int i; int prim; prims) 26 { 27 npts[i] = -len(primpoints(0, prim)); 28 } 29 prims = reorder(prims, argsort(npts)); 30 31 foreach(int i; int prim; prims) 32 { 33 int useNewColor = 1; 34 int pts[] = primpoints(0, prim); 35 for(int j = 0; j < ncolors; j++) 36 { 37 // See if we can use this color. 38 int useColor = 1; 39 foreach(int pt; pts) 40 { 41 if (mapping[j * npt + pt] != 0) 42 { 43 useColor = 0; 44 break; 45 } 46 } 47 48 if (useColor) 49 { 50 colors[i] = j; 51 foreach(int pt; pts) 52 { 53 mapping[j * npt + pt] = 1; 54 } 55 useNewColor = 0; 56 break; 57 } 58 } 59 if (useNewColor) 60 { 61 colors[i] = ncolors; 62 ncolors++; 63 resize(mapping, ncolors * npt); 64 foreach(int pt; pts) 65 { 66 mapping[colors[i] * npt + pt] = 1; 67 } 68 } 69 } 70 71 // Set prim color attribute. Colored prims will be skipped 72 // by the OpenCL pass. 73 foreach(int i; int prim; prims) 74 setprimattrib(geoself(), chs("colorname"), prim, colors[i]);