毕设论文 流渗模型的计算机模拟

图2－6 Al膜的随机电阻模拟[6]

6所示，可以将Al膜用随机电阻模拟，每个晶格用一个电阻代表，当电流通过时，只要存在一个通路就会有电流通过，可以将金属定向移动后留下空缺的过程认为是电阻的失效，随着电迁移的增加，电阻不断失效，当到一定程度后，突然会无电流通过，这时可以认为是Al线失效[6]。

图、PV图和LV图定义

图简介

许多科学以及工程问题可以转化为二维、三维的空间分割，如供货点供货区域划分问题以及无线通讯基站服务区域划分等问题。通讯基站或供货点位置可用孤立点表示，其空间划分多采用Voronoi模型，如图1所示,V图是计算几何中最重要的图维V图的定义：在N维空间R N上有n个点组成的点集S={p,p,p,…,p

图2-1 V图模型

2.3.2 LV图简介

为了更好的模拟实际的情况,对于模型进行了大量的改进,即加入了权重的概念。Laguerre-V oronoi图（LV图）就是一种施加权重的V图。LV图具体定义：设在N 维空间上有由n个球组成的集合G，G＝{c1,c2,…,c n},设r i,p i=（k1i,k2i, …,k N i）,分别是球c i的半径和球心坐标,定义空间一点p到球c i的距离d L(p,c i)为d L(p,c i)2=d(p,p i)2–r i2 ,于是可将满足R L（G,c i）＝{P∈R N| d(p,p i)< d(p,p j),|j≠i}的空间区域R L（G,c2）称为球c i的Laguerre-V oronoi区域(LV区域), 这样，N维空间被划分为N个区域和R L（G,c1）,R L（G,c2）,…,R L（G,c n）相应的边界，如此构成了LV图，或称power 图。此外，当区域R L（G,c i）和R L（G,c ij）有公共边时，可把两球的球心p i和p j连

图2-2 二维下的LV和LD图

中，给出了胞数为4096个胞的LV图形。

图2－3 实验所用的LV图（胞数为4096）

3.2 基于树的“连接/查找”的粗糙计算结果

本文采用树的“连接/查找”算法来进行程序的实现，分别在2D和3D上分别进行了仿真，得到在不同规模下的逾渗值，本试验中，每次进行完一次计算，当发生逾渗时，就将程序停止，并将结果记录到一个数组中保存起来，这样对于二维模型可以节省大约百分之五十的时间，对于三维模型效率就更加可观，在计算发生逾渗时增加一个break语句，可以省去发生逾渗后的循环过程，节约百分之八十的时间，这是非常可观的。

当然这个结果是比较粗糙的，本文在不同尺度下分别进行了一万次的计算，得到大量的样本值，这只是在微规则（microcanonical）下的结果(关于微规则与规则的概念以及精确计算将放在后面章节进行具体描述),本文计算这些这些逾渗值的平均值以及方差，其数值如表3-1所示：

表3－1 不同2D规模下的逾渗值p的平均值和方差

N 平均值方差

1024 0．5010 0．0414

2048 0．5008 0．0324

4096 0．4996 0．0253

8192 0．4996 0．0206

16384 0．4999 0．0161

从表中可以看出，在不同的规模下，逾渗值总是在0.5附近，这与pv图的逾渗值相差不大（本文同样可以计算pv图的逾渗值，在做LV图时只要让x，y坐标随机产生满足泊松分布即可，并让让r＝0，这时LV图将退化成为v图）。而方差随着规模

图3－4 2D下不同规模的逾渗值

的增大呈现减小的趋势，这也很好理解，在较小规模下p值分布更加不均匀可能添加

图3－5 3 D下的逾渗值分布

图3－6 2维下PV图不同规模的逾渗结果

3.4 性能分析

判断一个计算机程序好坏的一个重要的标准就是时间空间复杂度，通俗的说就是当一个程序在计算规模增加的情况下，运行的时间和数据增加的越少则复杂度越低，举例来说，如果计算时间如果与规模成幂或者指数增长时，那么这个算法的复杂度就较高，相应的对于工程计算来说，成本太高，不但不利于计算，也

图3－7 2D逾渗值计算的时间复杂度

图3－8 3D逾渗值计算的时间复杂度

图4－1 R.ziff作出的4种Rl（p）图形[8]

图4－2 一次试验的Rl（p）图形

在图4－3中，本文给出了胞数是1024时的Rl（n）的图形，可以看到这时图形还不够光滑，不仅存在突出的部分而且在峰顶附近有个很高的点。

图4－3 胞数是1024下的Rl（n）图形

在图4－4中给出了Rl（p）的图形，可以看出经过转换之后得到了一个更加平滑的图形，利用该图形的特性，峰值所对应的横坐标就是该规模下的逾渗值，本文利用matlab自带的max函数就可以找到对应的最大值的横坐标的位置。经计算这个规模下的逾渗值是0.5004(5)。

需要指出的是，这种转换并非是容易做到的，这里存在一个问题，看式（2）可知，在p较小且n较大时总有p n或者（1－p）N－n会非常非常小，以致于超过了matlab中的数值的极限，所有的编程软件中取值都有一个范围，matlab中当数据溢出时，比如当数非常小时，matlab认为这个数是0，而大于上限时则认为这个数是空，这给本文的运算带来了麻烦。

本文采用的解决方式是将最大的B（N,n,p）设为一个定值1，而二项分布在n ＝N×p时取得最大值，所以，将此时的B（N,n,p）设为1，所以由不同n下相互之

图4－4 在胞数是1024下的Rl（p）图形

图4－5 不同规模下的Rl(p)图形

下面本文给出3D的Rl（p）图形，二者很相似，但在原理上有一点区别：就是2D是记录的当x方向第一次发生逾渗时的座的占据百分数以及两个方向都发生逾渗时的座占据概率分别存入两个数组中，并以此来画出Rl（n）和Rl（p）的图形，而3D下是基于本文记录了在x，y轴向上均发生逾渗时的座占据概率以及在三个方向上都发生逾渗的概率来画出Rl（p）的图形。这事实上是两种不同的定义Rl的方式决定的。