最近在论文中看到有学者用改进的引力搜索算法解优化问题，有一个较好的效果，于是去了解了一下这个算法。
 
引力搜索算法(Gravitational Search Algorithm，GSA)是Esmat Rashedi等人在2009年提出的一种随机性启发式搜索算法，这种算法的灵感来自于牛顿的万有引力定律与运动定律：1.任意两个质点有通过连心线方向上的力相互吸引，该引力大小与它们质量的乘积成正比与它们距离的平方成反比。2.力使物体获得加速度。
 
在GSA中，质点被抽象成解空间的一个解，解之间存在一个相互吸引力，这个吸引力由解的质量与两个解之间的距离确定，质点的质量被抽象成解的评估函数值。在解空间中，每个解由其他解对其的吸引力获得加速度，质量更大(评估函数值更优)所提供的加速度更大，从而使解向更优解的方向移动。
 
看到这里，有些小伙伴会以为GSA是个与PSO差不多的算法，是的，GSA与PSO的外层框架是一致的，但是，它们最关键的粒子移动策略却是不同的：
 
1.在PSO中，粒子的移动只使用两个最佳位置pbest和gbest。但是在GSA中，智能体的方向是根据所有其他智能体或部分较优的智能体获得的总力来计算的。
 2.PSO使用一种内存来更新速度（由于pbest和gbest）。然而，GSA是无内存的，只有智能体的当前位置在更新过程中起作用。
 3.在PSO中，更新不考虑解之间的距离，而在GSA中，力与解之间的距离成反比。
 4.PSO模拟了鸟类的社会行为，而GSA的灵感来自于物理现象。
 
GSA的主要过程如下：
 
1. 确定搜索空间。
 2. 随机初始化个体种群。
 3. 对种群进行适应度评价。
 4. 更新引力常量G，更新种群中最好的个体 best与最差的个体worst，更新个体的质量。
 5. 计算每个个体在不同方向上的总引力。
 6. 计算个体的加速度，基于此更新个体的速度。
 7. 根据速度更新个体在解空间中的位置。
 8. 重复步骤3-7直到达到停止标准。
 
算法流程：
 
 
具体的个体属性更新公式与引力计算公式建议参考 GSA: A Gravitational Search Algorithm这篇文章。
 
例： min z = x1 ** 2 + x2 ** 2
 s.t. -10 <= x1 <= 10 , -10 <= x2 <= 10
 
# Gravitational Search Algorithm
import numpy as np
import random as rd
from math import exp, sqrt

'''min z = X1 ** 2 + X2 ** 2
   s.t.   -10 <= x1 <= 10
          -10 <= X2 <= 10   '''

def init(n):
    position, velocity = [], []
    for i in range(n):
        X1 = rd.uniform(-10, 10)
        X2 = rd.uniform(-10, 10)
        V1 = rd.uniform(-3, 3)
        V2 = rd.uniform(-3, 3)
        position.append([X1, X2])
        velocity.append([V1, V2])
    return position, velocity

def objFuntion(x1, x2):
    return x1 ** 2 + x2 ** 2

def fitnessEva(position):
    fitness = []
    for i in range(len(position)):
        fitness.append(objFuntion(position[i][0], position[i][1]))
    return fitness

def findBestAndWorst(position):
    return min(fitnessEva(position)), max(fitnessEva(position))

def calculateMass(fitness):
    mass = []
    Mass = []
    for i in range(len(fitness)):
        mass.append((fitness[i] - max(fitness)) / (min(fitness) - max(fitness)))
    for i in range(len(mass)):
        Mass.append(mass[i] / sum(mass))
    return Mass

def calculateAcceleration(position, Mass, G, topK):
    acceleration = []
    Fi0, Fi1 = 0, 0

    for i in range(len(position)):
        for j in range(len(position)):
            if i != j and j in topK:
                Fi0 += rd.random() * G * ((Mass[i] * Mass[j]) / (calculateDistance(position[i], position[j]) + r)) * (
                            position[j][0] - position[i][0])
                Fi1 += rd.random() * G * ((Mass[i] * Mass[j]) / (calculateDistance(position[i], position[j]) + r)) * (
                            position[j][1] - position[i][1])
        if Mass[i] != 0:
            acceleration.append([Fi0 / Mass[i] / 10, Fi1 / Mass[i] / 10])   #这里除10是为了避免粒子的加速度过大
        else:
            acceleration.append([10, 10])
        Fi0 = 0
        Fi1 = 0
    return acceleration

def findTopK(fitness, K):
    topK = []
    dic = {}
    for i in range(len(fitness)):
        dic[i] = fitness[i]
    fitness.sort()
    for i in range(K):
        topK.append(list(dic.keys())[list(dic.values()).index(fitness[i])])
    return topK


def updateVelocityAndPosition(acceleration, position, velocity):
    for i in range(len(velocity)):
        velocity[i][0] = rd.random() * velocity[i][0] + acceleration[i][0]
        velocity[i][1] = rd.random() * velocity[i][1] + acceleration[i][1]

        position[i][0] = position[i][0] + velocity[i][0]
        position[i][1] = position[i][1] + velocity[i][1]

def calculateDistance(p1,p2):
    return sqrt((p1[0] - p2[0]) ** 2 + (p1[1] - p2[1]) ** 2)

def checkPosition(position):
    for i in range(len(position)):
        if position[i][0] < -10:
            position[i][0] = -10
        elif position[i][0] > 10:
            position[i][0] = 10
        if position[i][1] < -10:
            position[i][1] = -10
        elif position[i][1] > 10:
            position[i][1] = 10

if __name__ == '__main__':
    G = 100
    r = 1
    K = 50

    iterx, maxIterx = 0, 50
    position, velocity = init(50)
    while iterx < maxIterx:
        fitness = fitnessEva(position)               #适应性评估
        G = G * exp(-20 * iterx / maxIterx)          #更新引力常量
        Mass = calculateMass(fitness)                #更新粒子质量
        topK = findTopK(fitness, K)                  #找出适应度更优的前K个粒子
        acceleration = calculateAcceleration(position, Mass, G, topK)          #计算粒子加速度
        updateVelocityAndPosition(acceleration, position, velocity)            #根据加速度更新速度与位置
        checkPosition(position)                                                #检查粒子是否冲出了解空间
        iterx += 1
        K = K - iterx                                #更新K值
    print(min(fitnessEva(position)))
    print(position[fitnessEva(position).index(min(fitnessEva(position)))])
 
参考文献：
 Rashedi E, Nezamabadi-pour H, Saryazdi S. GSA:a gravitational search algorithm[J].Information Science, 2009, 179 (13) :2232-2248.

 参数说明：
 $M:地球总质量$
 $dm:地球上的质量微元$
 $m^{\prime }:地球外部空间上的一个质点$
 $p:m与dm之间的距离$
 万有引力大小为：
 $f=\frac{Gdm\ast m^{\prime }}{p^2}$
 
引力位函数
 
$若外有引力做工(从无穷远移动到半径p处)为:A={\int }_{\infty }^p-\frac{Gdm\ast m}{p^2}=\frac{Gdm\ast m^{\prime }}{p}$
 根据能量守恒，引力做工必然等于$m^{\prime }$减少的位能（势能），若单位质量的点(从无穷远移动到半径p处)，势能必然减少，将减少量定义为位函数:
 $$dV=\frac{Gdm}{p}$$
 
地球外部空间的调和函数
 
地球产生的各质量微元函数之和为：
 $$V={\int }_{M}dV=G{\int }_M\frac{dm}{p}（公式1）$$
 又因为：
 $$▽V=\frac{{\partial }^{2}V}{\partial x^2}+\frac{{\partial }^{2}V}{\partial y^2}+\frac{{\partial }^{2}V}{\partial z^2}=0(证略)$$
 所以$V$是在地球质量M的外部空间上是调和函数
 
勒让德多项式的生成函数
 
 
参数说明：
 $m^{\prime }：的球坐标为(r,\theta ,\lambda )$
 $dm：的球坐标为(r,{\theta }^{\prime },{\lambda }^{\prime })$
 $向量R与r夹角为\phi$
 $m{0}^{\prime }是m^{\prime }在球面上的投影$
 在,$o,dm,m^{\prime }$构成的三角形中，根据余弦定理,可以得到：
 $$p^2=r^2+R^2-2Rrcos\phi =r^2(1-2ax+a^2)(a=\frac{R}{r},x=cos\phi )$$
 再开平方，取倒数得到：
 $$\frac{1}{p}=\frac{1}{r}(1-2ax+a^2{)}^{-\frac{1}{2}}（公式2）$$
 又因为：
 $$(1-2ax+a^2{)}^{-\frac{1}{2}}=\sum _{n=0}^{\infty }{P}_n(x)a^n(勒让德多项式的生成函数)（公式3）$$
 所以$(1-2ax+a^2{)}^{-\frac{1}{2}}$为勒让德多项式的生成函数
 
球函数的加法公式
 
 
根据球面的三角余弦定理有:$cos\phi =cos\theta cos{\theta }^{\prime }+sin\theta sin{\theta }^{\prime }cos(\lambda -{\lambda }^{\prime })$
 带入$$P_n(cos\phi )=\sum _{k=0}^n\frac{2(n-k)!}{(1+{\delta }_k)(n+k)!}{P}_n^k(cos\theta )P_n^k(cos{\theta }^{\prime })cosk{\lambda }^{\prime }+P_n^{k}cos\theta sin(k\lambda )P_n^k(cos{\theta }^{\prime })(sink{\lambda }^{\prime })(公式4)$$
 
利用德让勒多项式求解V（调和函数）
 
将勒让德多项式的生成函数（公式(3)），将带入公式（2），再将公式2带入地球外部空间的调和函数公式（1）。
 得到:
 $R_e:旋转椭圆球的长半径$
 $u=GM:地球引力常数$
 ${\delta }_k=\{\begin{array}{ll}1& k=0\\ 0& k\ne 1\end{array}$
 $$V=\frac{G}{r}\sum _{n=0}^{\infty }(\frac{R_e}{r}{)}^n{\int }_M(\frac{R}{R_e}{)}^n{P}_n(cos\phi )dm$$
 将公式（4）带入:
 得到:
 $$V=\frac{u}{r}\sum _{n=0}^{\infty }(\frac{R_e}{r}{)}^n\sum _{k=0}^n(C_n^{k}cosk\lambda +S_n^{k}sink\lambda )P_n^k(cos\theta )$$
 $$C_n^k=\frac{2(n-k)!}{M(1+{\delta }_k)(n+k)!}{\int }_M(\frac{R}{R_e}{)}^n{P}_n^k(cos{\theta }^{\prime })cosk{\lambda }^{\prime }dm$$
 $$S_n^k=\frac{2(n-k)!}{M(1+{\delta }_k)(n+k)!}{\int }_M(\frac{R}{R_e}{)}^n{P}_n^k(cos{\theta }^{\prime })sink{\lambda }^{\prime }dm$$
 已知直角坐标系到球体坐标系的转换关系如下：
 $$\{\begin{array}{l}x=Rsin{\theta }^{\prime }cos{\lambda }^{\prime }\\ y=Rsin{\theta }^{\prime }sin{\lambda }^{\prime }\\ z=Rcos{\theta }^{\prime }\end{array}$$
 刚体的张量定义如下:刚体的惯性张量及其物理意义
 $$I=\left[\begin{array}{ccc}{I}_x& -I_{xy}& -I_{xz}\\ -I_{xy}& I_y& -I_{yz}\\ -I_{xz}& -I_{yz}& I_{xz}\end{array}\right]={\int }_M\left[\begin{array}{ccc}{y}^2+z^2& -xy& -xz\\ -xy& x^2+z^2& -yz\\ -xz& -yz& x^2+y^2\end{array}\right]dm$$
 可以得到：
 $C_0^0=1$
 $C_0^1=\frac{1}{R_e}(\frac{1}{M}{\int }_{M}zdm)$
 $C_1^1=\frac{1}{R_e}(\frac{1}{M}{\int }_{M}xdm)$
 $S_1^1=\frac{1}{R_e}(\frac{1}{M}{\int }_{M}ydm)$
 $C_2^0=\frac{1}{M{R}_e^2}({\int }_M{z}^2-\frac{x^2+y^2}{2}dm)=\frac{1}{M{R}_e^2}(\frac{I_x+I_y}{2}-I_z)$
 $C_2^1=\frac{1}{M{R}_e^2}{\int }_{M}xzdm=\frac{I_{xz}}{M{R}_e^2}$
 $C_2^2=\frac{1}{4M{R}_e^2}({\int }_M{x}^2-y^{2}dm)=\frac{I_y-I_x}{4M{R}_e^2}$
 $S_2^1=\frac{1}{M{R}_e^2}({\int }_{M}yzdm)=\frac{I_{yz}}{M{R}_e^2}$
 $S_2^2=\frac{1}{4M{R}_e^2}({\int }_{M}xydm)=\frac{I_{xy}}{2M{R}_e^2}$
 如果定义直角坐标与地球惯性惯性主轴重合，则：$I_{xy}=I_{yz}=I_{xz}=0$
 
对于实际应用时，地球坐标一边选择坐标原点与地球质心重合，oz与地球自转平行，ox轴在赤道面上且指向0度经线，这时坐标轴往往与惯性主轴不重合。
 利用三角函数恒等式：
 得到:
 
     
      
       
        
         V
        
        
         =
        
        
         
          u
         
         
          r
         
        
        
         {
        
        
         1
        
        
         −
        
        
         
          ∑
         
         
          
           n
          
          
           =
          
          
           1
          
         
         
          ∞
         
        
        
         (
        
        
         
          
           R
          
          
           e
          
         
         
          r
         
        
        
         
          )
         
         
          n
         
        
        
         [
        
        
         
          J
         
         
          n
         
        
        
         
          P
         
         
          n
         
        
        
         (
        
        
         c
        
        
         o
        
        
         s
        
        
         θ
        
        
         )
        
        
         +
        
        
         
          ∑
         
         
          
           k
          
          
           =
          
          
           1
          
         
         
          n
         
        
        
         
          J
         
         
          n
         
         
          k
         
        
        
         
          P
         
         
          n
         
         
          k
         
        
        
         (
        
        
         c
        
        
         o
        
        
         s
        
        
         θ
        
        
         )
        
        
         c
        
        
         o
        
        
         s
        
        
         k
        
        
         (
        
        
         λ
        
        
         +
        
        
         
          λ
         
         
          n
         
         
          k
         
        
        
         )
        
        
         ]
        
        
         }
        
       
       
        V=\frac{u}{r}\{1-\sum_{n=1}^\infty(\frac{R_e}{r})^n[J_nP_n(cos\theta)+\sum_{k=1}^nJ_n^kP_n^k(cos\theta)cosk(\lambda+\lambda_n^k)]\}
       
      
     V=ru​{1−n=1∑∞​(rRe​​)n[Jn​Pn​(cosθ)+k=1∑n​Jnk​Pnk​(cosθ)cosk(λ+λnk​)]}
 $$J_n=-C_n^0,J_n^k=\sqrt{(C_n^k{)}^2+(S_n^k{)}^2},{\lambda }_n^k=-arctan(\frac{S_n^k}{C_n^k})/k$$
 $J_2=C_2^0=\frac{1}{M{R}_e^2}(\frac{I_x+I_y}{2}-I_z)(动力扁率：反映赤道与极轴上转动惯量的差别)$
 $\frac{u}{r}:球形地球引起的引力位$
 在实际应用中因为其他系数比$J_2$小三个数量级，所以只需考虑$\frac{u}{r}和J_2$的影响，可以将$V$近似看为：
 $$V=\frac{u}{r}[1-\frac{J_2{R}_e^2}{2r^2}(3co{s}^2\theta -1)]（只考虑主谐项）$$
 
卫星在惯性坐标系下的运动方程（只考虑主谐项）
 
已知
 $位移矢量：r={\left[\begin{array}{ccc}x& y& z\end{array}\right]}^T$
 $引力：f={\left[\begin{array}{ccc}{f}_x& f_y& f_z\end{array}\right]}^T$
 $cos\theta =\frac{z}{r}$
 对x,y,z求2阶片导，求得加速度方程为：
 $$\ddot{r}=-\frac{u}{r^3}[I+\frac{3}{2}{J}_2(\frac{R_e}{r}{)}^2(D-5(u_r\ast u_p)I)]r$$
 其中:$D=diag\left(\begin{array}{ccc}1& 1& 3\end{array}\right)$
 $u_r表示r上的单位矢量$
 $u_p表示自转轴上的的单位矢量$
 $u_r\ast u_p表示u_p与u_r之间夹角的余弦值$

​公众号原文链接，欢迎关注基于引力模型的城市引力分析-城市经济联系图城市引力模型—制作城市经济联系图https://mp.weixin.qq.com/s/EuLjbC8UbrRZh7dJHrVNvg
 
还是之前写的城市联系度咋实现，我稍微改了一点点，最后用Qgis调整了一下图片。
 
 
 
一、瞎扯
 
在知网上一搜“基于引力模型的空间格局研究”，有各省份相关的研究，也是通过GDP、人口等一切你能自圆其说的社会经济数据，利用ArcGIS来做各城市间的联系度。如：
 
 
（图片来自《中原经济区城市间相互作用时空格局演变研究》）
 
无论是用什么样的数据，原理都是基牛顿力学的万有引力公式为基础。
 
 
Tinbergen（1962）和Poyhonen（1963）对其在经济学领域做了发展、延伸，提出了一个比较完整且简便的经济学模型—引力模型。这个模型认为两个经济体之间的单项贸易流量与它们各自的经济规模（一般用GDP来表示）成正比，与它们之间的距离成反比。
 
两个城市之间的引力模型公式通常可表现为：
 
 
 
其中，F为两个城市之间的引力，Q为两个城市的规模，d为两个城市之间的距离，a为常量，不影响，可取值为1。
 
 
城市规模指标Q可以是经济规模、人口规模、消费规模、进出口规模等等单项或综合；距离d可以采用实际的公路里程数，也可以采用两城市之间的空间直线距离（欧氏距离）等单项或综合。
 
可以采取城市综合指标处理的方式，一般采用的是求几何平均。例如选取4个指标的几何平均值作为城市规模指标。
 
 
在GIS里就是把数据弄进属性表里，然后用字段计算器用引力公式计算结果，最后分级表现。
 
 
要实现城市间的经济联系图，需要将每个城市分别作为起点和终点两两组合得到城市之间的连线，将每条线都附上城市之间引力值（通过公式求出的值）便能得到城市联系网络图。
 
那么要在ArcGIS里实现，需要做的事有两件。
 
 
每个城市的x,y坐标点
 
 
城市之间两两组合排列
 
之后使用xy转线工具，得到城市之间的联线。最后将计算出的引力值，将其连接至线上即可。
 
 
二、操作步骤
 
以2019年山东省16个地级市为例，选取GDP作为城市规模指标（单一，肯定不是很严谨，仅供参考）。
 
 
Q为GDP，d为距离，a为常量，不影响，可取值为1。
 
将山东省16地市面要素转为点要素，在属性表里添加了xy字段，通过计算几何计算得出xy坐标。各市GDP数据连接进来。然后是数据的投影，不然计算距离用地理坐标系是不行的。
 
 
用【点距离】工具来计算，注意这里输出表路径选择默认，我选择其他路径会报错，不知道什么原因。
 
 
打开输出的表，每个FID都代表一个城市，input_FID相当于起点城市，near_FID就相当于终点城市，distance便是两个城市间的距离。
 
上表中的FID与下面这个表中的input_FID，near_FID是对应的。需要把它们关联起来，不然看数字谁知道是哪个城市。
 
 
可以在GIS里进行连接和关联，也可以直接复制到excle表格里实现。我jiao的用excle简单些，所以我选择excle来做。
 
复制上面两步的两个属性表到excle（Ctrl+A全选-右键FID左边空白-复制所选项），在excle里把用不着的字段删掉就可以。
 
 
复制出来得到这两个表
 
 
在第一个表格里新增加6个表头，分别是or代表起点城市，ox、oy代表起点城市坐标；des代表终点城市，dx、dy代表终点城市坐标。
 
 
接下来用到的就是excle里的vlookup函数了。简单说一下vlookup的意思：
 
Lookup_value是我们需要对应的那一列，先做的起点城市，先选input_fid。
 
Table_array是对应的哪个表，选中sheet2全部。
 
Col_index就是对应的sheet2表里的哪一列，城市名称是第三列，输入3即可。
 
 
剩下的表头填充都是这方法，填充进去即可，不再演示。最后得到包含起点终点城市名，xy坐标，起终点GDP数据的表格。
 
 
接下来用开头的公式进行计算
 
 
最开始计算的距离单位是米，这里把距离除1000，换算成千米。
 
Value（引力值）=ogdp（起点城市GDP）*dgdp(终点城市GDP)/距离的平方
 
 
打开ArcGIS里XY转线工具，依次输入选项，ID选择计算出来的引力值VALUE。
 
 
这样我们初步的结果就出来了，可以进一步对其可视化。最后大致是这个样子
 
 
Qgis作图比Arcgis的功能好多了，这次就放到Qgis里面试了一下，操作与Arcgis都是相通的，加载一下在线底图和图层文件，同样对线进行符号化显示就可以。
 
 
可以看到，半岛地区，青岛与潍坊和烟台之间引力较大。西部就是济南与淄博，其次是德州和泰安。仅用了GDP这个单一指标进行试验，很难相对全面考察联系强度，初始重力模型太简单，也没去看很多论文再改进，如果仅用作平时项目分析之类的我觉得可以用用。
 
 
最后，Qgis加载在线底图的插件下载地址：
 
https://plugins.qgis.org/plugins/HCMGIS/
 
Arcgis还有一个OD插件工具可以把这种直线变为弧度线，我没试过，链接也放在这，有需自取：
 
https://pan.baidu.com/s/1UIuBIOWzVnPmrhXaZF7ggQ
 
提取码：c0es
 
 
END
 
如果在实验过程中有什么地方遇到问题可以加我，竭力解答。需要练习数据也请后台我，公主号【一点规划】。
 
如果可以的话，希望能够点个关注转发分享，点个在看并且点个赞~~
 

 
 
matlab绘制引力场
 
 您在一个街区。 有风，小雨正在下。 人们来回穿梭，携带着大量的粉笔和防水布。 许多人蹲在街上，他们的粉笔艺术受苦并诅咒着雨。 这条街向南延伸。 一条肮脏的混凝土楼梯井通向街道下方。 
  
 
 
 >下去 
 
 
 楼梯的底部是一扇敞开的门，通向一个黑暗的房间。 黑暗似乎从楼梯间底部的门口和游泳池溢出。 楼梯向上。 一个没有照明的门向东打开。 
 
 
 >往东走 
 
 
 您输入。 有了这么少的光线，您看不到多少东西。 您举起手，几乎无法辨认其形状。 当您抽出手时，远处墙壁上的一小束光突然爆炸，使您的视线陷入彩虹。 
 
 

 

  
 
 
 
 欢迎来到DiMo 
 
 
 Digital Motion（DiMo）是GeekSpark的年度互动艺术展览，该展览每年在北卡罗来纳州罗利的SparkCon期间举行。 
 
 
 从表面上看，这是您挥舞手臂与巨大的投影屏进行交互的事件（今年，玩家挥舞着发光的警棍与计算机生成的图形进行交互）。 从更深层次上讲，GeekSpark的目标是使访客充满乐趣，然后使他们能够探索其背后的技术。 为此，DiMo背后的所有代码都是开源的。 
 
 
 2014年是我作为志愿者的第一年，与所有人一起工作是一次非常有趣的教育。 以下是制作这些展览品时团队合作的简短描述： 
 
 
 红帽工程师Ian Hands创建了一个光传感器系统，使玩家可以通过挥动警棍来与显示器互动。 共有三个可视化，所有这些可视化都是围绕光传感器系统构建的。 显示屏从一次可视化显示每五分钟旋转一次。 
 
 
 乔纳森·里皮（Jonathan Rippy）创作了《 康威的人生游戏 》的互动版本。 
 
 

 

  
 
 
 
 Wray Bowling创建了一个有趣的三人游戏，结合了多力多滋，乡巴佬，几何和苏打水。 
 
 

 

  
 
 
 
 我创建了一个粒子物理sim卡，允许玩家与重力场（和漂亮的颜色）进行交互。 
 
 

 

  
 
 
 
 在探讨技术细节之前，请查看此现场演示 。 
 
 
 该演示使您可以使用鼠标而不是发光警棍来玩粒子模拟。 需要现代的网络浏览器，最好是Firefox或Chrome。 按空格键在播放器之间切换，然后使用设置面板滑块创建全新的可视化效果。 
 
 

 

  
 
 
 
 DiMo的起源 
 
 
 去年五月，伊恩（Ian）与我联系，为GeekSpark做贡献。 经过一连串的休息室讨论和具有象征意义的餐巾素描后，我开始研究粒子模拟，而他正在研究传感器系统。 Jonathan和Wray从Red Hat外部慷慨地自愿捐款。 
 
 
 这个想法是：将发光的警棍散发给玩家，然后他们可以挥手与酷炫的墙面显示器互动。 工作原理：一台在天花板附近装有广角网络摄像头的计算机将在每个播放器的位置上保持标签，并将这些位置传递给可视化效果，从而使播放器部件在屏幕上移动。 
 
 
 DiMo的三种可视化效果是由HTML，CSS和JavaScript组成的常规网页。 之所以选择将Web用于此项目，是出于以下几个原因： 
 
 
 较低的进入壁垒是一个很大的障碍。 只需几分钟和少量指导，即使没有编程经验的人也可以开始自己做。 图形编程在所有类型的编程中都是独一无二的。 通常，方程式的意外变化（甚至是错字！）都会带来非常酷的惊喜。 “仅需尝试”就能获得比其他任何编程分支更大的收益。 
 活动结束后共享页面就像单击链接一样简单。 与其他平台相比，网络更容易将您的像素摆在人们面前。 
 Web本质上是开放的。 每个（桌面）网络浏览器都有一个“查看源代码”按钮，因此任何有意愿的人都可以找到代码。 
 
 
 接下来，DiMo 2015！ 
 
 
 您想制作自己的DiMo显示器吗？ 保持联系！ 2015年DiMo活动的计划即将开始，因此没有比现在更好的参与时间了！ 对Raleigh而言是本地的，这很好，但是并不需要。 不过，如果您可以降低SparkCon的积分，您将获得加分。 看到所有陌生人都喜欢自己制作的东西，这无比喜悦。 
 
 

 

  
 
 
 
 链接和其他好奇心： 
 
 
 传感器系统（输入）和粒子模拟（渲染器）的源代码 
 粒子模拟的视频 
 设置晚会的间隔拍摄视频 
 演示如何使用控件创建新效果的视频 
 有关粒子模拟的更多信息 
 Wray Bowling的图片组合 
 
 

  
翻译自: https://opensource.com/life/15/2/sparkcon-geekspark-digital-motion-exhibit
 
 
matlab绘制引力场