（计算机软件与理论专业论文）多核环境下受限voronoi图的研究和实现.pdf

摘要 摘要 v o r o n o i 图作为空间离散划分的基本几何结构单元，是计算几何领域重要的 研究方向之一。v o r o n o i 图因其垂直平分的特性，特别适用于利用有限体积法求 解如液体间流动、热传导等具有某种守恒性质的物理问题，但是在这些应用中， 对v o r o n o i 图的生成有特殊的限定条件，以油藏模拟为例，生成的v o r o n o i 网格 单元不能跨越边界、断层等物理障碍，对于一些特殊的点，如水平井点和垂直井 点必须为v o r o n o i 网格的中心点。由于实际应用中对受限v o r o n o i 图的需求，以 及受限情况下构造高质量的v o r o n o i 网格需要进行大规模数据的计算处理，因此 有必要对受限v o r o n o i 网格的构造及其并行化算法进行研究和实现。 在现有的单机单核的计算能力下，处理大规模数量的网格生成，算法效率低 下，计算时间无法接受，而在大型并行处理机进行数据处理，成本高，用户操作 不便。随着硬件处理器的快速发展，以异构多核为特征的众核处理器已成为目前 高性能计算的主流趋势，因此，结合实际需求，本文提出了在多核环境下利用 m p i 并行编程语言实现高效、可用的网格生成并行算法，具体研究成果和内容包 括以下几个方面： ( 1 ) 根据d e l a u n a y 图是v o r o n o i 图的对偶图的性质，使用间接法实现 v o r o n o i 图构造算法，生成网格结构饱满，质量高，但运行速度较慢。 ( 2 ) 基于算法并行化考虑，采用最优算法一分治法来实现带有限定情况的 受限v o r o n o i 图的生成。根据两个限定条件( 限定点和限定线原则) ， 对各种限定情况做了分析处理，为后期生成高质量的网格进行合理布 点。 ( 3 ) 利用分治法构造大规模点集的v o r o n o i 网格，同时设计修改了数据结 构双向循环链表( d o u b l y - c o n n e c t e de d g el i s t ，d c e l ) ，使之符合并行 处理和限定v o r o n o i 复杂边界的特点。 ( 4 ) 使用m p i 和c 语言，采取逐步合并的整体并行方案，对合理布点的 点集生成和v o r o n o i 构造做了并行化处理，在大型惠普并行机和异构 多核单机上分别进行了运算，得到实验结果。 ( 5 ) v o r o n o i 对偶图d e l a u n a y 的转换实现，最终对数据进行处理实现图形 的实时显示。 摘要 ( 6 ) 最后，本文对基于v o r o n o i 图技术未来的研究方向，进行了展望。 本文的研究成果可应用于模拟流体流动，建模绘制等领域，结合多核处理器 和并行计算高速发展的时代特点，用于各种模拟问题的实时计算。 关键词：受限v o r o n o i 、d e l a u n a y 三角剖分、分治算法、并行计算、m p i 、 异构多核 a b s t r a c t a b s t r a c t a sab a s i cd a t as t r u c t u r ee l e m e n to fd i s c r e t es p a c ed i v i s i o n ，v o r o n o im e s hi sa n i m p o r t a n tr e s e a r c hd i r e c t i o ni nt h ef i e l do fc o m p u t a t i o n a lg e o m e t r y b e c a u s eo ft h e c h a r a c t e r i s t i co fp e r p e n d i c u l a rb i s e c t o r , v o r o n o im e s he s p e c i a l l ys u i t a b l ef o rs o l v i n g s o m ep h y s i c a lp r o b l e m sr e l a t e dt oc o n s e r v a t i o nn a t u r eb yt h em e t h o do ff i n i t ev o l u m e ， s u c ha sf l u i df l o w , h e a tc o n d u c t i o na n ds oo n h o w e v e r , i nt h e s ea p p l i c a t i o na r e a s ，t h e g e n e r a t i o no ft h ev o r o n o id i a g r a mh a sas p e c i a lq u a l i f i c a t i o n ，t a k et h ee x a m p l eo f r e s e r v o i rs i m u l a t i o n ，t h ev o r o n o im e s hc a nn o tc r o s ss o m ep h y s i c a lb a r r i e r s ，l i k et h e b o r d e ra n df a u l t s f o rs o m ep a r t i c u l a rp o i n t s ，s u c ha sh o r i z o n t a lw e l la n dv e r t i c a lw e l l m u s tb et h ec e n t e ro ft h ev o r o n o im e s h c o n s i d e rt h er e q u i r e m e n to fc o n s t r u c t i n g l i m i t e dv o r o n o im e s hi ns o m ep r a c t i c a la p p l i c a t i o n sa n dt h ep r o b l e mo fc o n s t r u c t i n g h i g h q u a l i t y v o r o n o im e s hu n d e rt h ec o m p l e xc o n s t r a i n e dc o n d i t i o n s t h r o u g h l a r g e s c a l ed a t ap r o c e s s i n g ，r e s e a r c ha n di m p l e m e n t a t i o no fl i m i t e dv o r o n o im e s h g e n e r a t i o na n di t sp a r a l l e la l g o r i t h m si sn e c e s s a r y u n d e rt h ec o m p u t i n ga b i l i t yo fs i n g l e - c o r e ，u s i n gl a r g e - s c a l ed a t at og e n e r a t e l i m i t e dv o r o n o im e s hi sl o we f f i c i e n c ya n di t sc a l c u l a t i o nt i m ec a nn o tb ea c c e p t e d ， h o w e v e r , u s i n gt h el a r g e - s c a l ep a r a l l e lp r o c e s s o rf o rd a t ap r o c e s s i n gn e e d sh i g hc o s t ， a n di t sn o tc o n v e n i e n tt oo p e r a t i o n w i t ht h er a p i dd e v e l o p m e n to fh a r d w a r ea n d p r o c e s s o r ，h e t e r o g e n e o u sm u l t i c o r eh a sb e c o m et h em a i nt r e n do fh i g hp e r f o r m a n c e c o m p u t i n g i nt h i sp a p e r , w ep r e s e n taa v a i l a b l ea n de f f i c i e n t v o r o n o im e s h g e n e r a t i o np a r a l l e la l g o r i t h mb yu s i n gm p ip a r a l l e lp r o g r a m m i n gl a n g u a g eu n d e rt h e e n v i r o n m e n to fm u l t i - c o r e ，t h es p e c i f i cr e s e a r c ha n dc o n t e n ta st h ef o l l o w s ： ( 1 ) i m p l e m e n t a t i o nt h ea l g o r i t h mo fc o n s t r u c t i n gv o r o n o im e s hb yt h ei n d i r e c t m e t h o d ，w h i c hi sh i g h q u a l i t yb u tl o we f f i c i e n c y ( 2 ) c o n s i d e rt h ep a r a l l e la l g o r i t h m ，w eu s et h eo p t i m a la l g o r i t h mo fd i v i d ea n d c o n q u e rt o a c h i e v el i m i t e dv o r o n o im e s hg e n e r a t i o n a c c o r d i n gt ot h et w o q u a l i f i c a t i o n s ( t h ep r i n c i p l e so fl i m i t e dp o i n ta n dl i m i t e dl i n e ) ，w em a k e a n a l y s i sa n dp r o c e s s i n g o ne a c hl i m i t e ds i t u a t i o n ，w h i c hi s p r e p a r e f o r c o n s t r u c t i n gh i g h - q u a l i t ym e s hg e n e r a t i o nl a t t e r i i i a b s t r a c t ( 3 ) c o n s t r u c t i n gl i m i t e dv o r o n o im e s ho fl a r g e - s c a l ep o i n ts e tb yo p t i m a l a l g o r i t h m d i v i d ea n dc o n q u e ra l r g o r i t h m ，w ea l s od e s i g na n dm o d i f yt h e p a r a l l e ld a t as t r u c t u r e d o u b l y c o n n e c t e de d g el i s t ，w h i c hi se a s yt oc o n f o r m p a r a l l e lp r o c e s s i n ga n dc o m p l e xf e a t u r e so fc o m p l e xl i m i t e db o u n d a r y ( 4 ) u s i n gm p ia n dcp r o g r a m m i n gl a n g u a g e ，a c c o r d i n gt ot h eg r a d u a l l ym e r g ei n t h ew h o l ep a r a l l e l p r o g r a m ，w ed ot h ep a r a l l e lp r o c e s s i n go fr e a s o n a b l e d i s t r i b u t i o no fp o i n ta n dc o n s t r u c t i n gl i m i t e dv o r o n o im e s hg e n e r a t i o n w e r u nt h ep r o g r a mu n d e rt h ee n v i r o n m e n to fo nt h eh pl a r g e s c a l ep a r a l l e l m a c h i n ea n dh e t e r o g e n e o u sm u l t i c o r em a c h i n e ，a n dg e ti t se x p e r i m e n t a l r e s u l t s ( 5 )a c h i e v et h ev o r o n o id i a g r a mt r a n s f o r mt oi t sd u a ld i a g r a m - d e l a u n a yd i a g r a m a n dr e a l i z er e a l t i m eg r a p h i c sd i s p l a yb yu l t i m a t ed a t ap r o c e s s i n g ( 6 ) f i n a l l y , t h ef u t u r er e s e a r c hd i r e c t i o no fv o r o n o id i a g r a mi sp r e s e n t e d b a s e do nt h ec o n t e n to ft h i sp a p e r , w ec a na p p l yo u rr e s e a r c hr e s u l t so nt h ef i e l d s o ff l u i df l o ws i m u l a t i o n ，m o d e l i n gd r a w i n ga n ds oo n c o m b i n e dw i t hm u l t i c o r e p r o c e s s o r sa n dt h er a p i dd e v e l o p m e n to fp a r a l l e lc o m p u t i n g ，w ec a nu s et h i sp a r a l l e l a l g o r i t h mt os o l v em a n yp r o b l e m so fr e a l t i m ec o m p u t i n gr e l a t e dt os i m u l a t i o n k e yw o r d s ：l i m i t e dv o r o n o i ，d e l a u n a yt r i a n g u l a t i o n ，d i v i d ea n dc o n q u e r a l g o r i t h m ，p a r a l l e lc o m p u t i n g ，m p il a n g u a g e ，h e t e r o g e n e o u s m u l t i c o r e i v 中国科学技术大学学位论文原创性声明 本人声明所呈交的学位论文，是本人在导师指导下进行研究工作所取得的成 果。除己特别加以标注和致谢的地方外，论文中不包含任何他人已经发表或撰写 过的研究成果。与我一同工作的同志对本研究所做的贡献均已在论文中作了明确 的说明。 作者签名： 中国科学技术大学学位论文授权使用声明 作为申请学位的条件之一，学位论文著作权拥有者授权中国科学技术大学拥 有学位论文的部分使用权，即：学校有权按有关规定向国家有关部门或机构送交 论文的复印件和电子版，允许论文被查阅和借阅，可以将学位论文编入有关数据 库进行检索，可以采用影印、缩印或扫描等复制手段保存、汇编学位论文。本人 提交的电子文档的内容和纸质论文的内容相一致。 保密的学位论文在解密后也遵守此规定。 、 口公开吻保密( 三年) 作者签名：羔盛 签字日期： 逸f ! ：圭! 导师签名：弛 第1 章绪论 第一章绪论 计算几何( c o m p u t a t i o n a lg e o m e t r y ) 作为计算机科学的一个分支，主要研究 解决几何问题的算法，其研究内容涉及到凸包( c o n v e x h u l l ) 、v o r o n o i 图( v o r o n o i d i a g r a m s ) 、三角化( t r i a n g u l a t i o n ) 、直线求交( a l l - i i n e - i n t e r s e c t i o n s ) 、多边形划分 ( p o l y g o nd i v i s i o n ) 、区域查找( r e g i o n a ls e a r c h ) 、距离计算( d i s t a n c ec a l c u l a t i o n ) 、 可见性计算( v i s i b i l i t y ) 、路径规划( p a t hp l a n n i n g ) 、点定位( p o i n tp o s i t i o n i n g ) 等内 容。在现代工程和数学领域，计算几何在图形学、机器人技术、超大规模集成电 路设计和统计等诸多领域有着十分重要的应用。 v o r o n o i 图作为计算几何领域重要的几何基元，是计算几何的基本几何结构 单元，因其广泛的应用领域和重要的理论价值引起了众多学者的研究兴趣，发表 相关v o r o n o i 研究领域的文章占计算几何领域的1 5 左右。简而言之，v o r o n o i 图 一直是计算机几何学科中的核心问题之一，也是当前国内外计算机几何研究的热 点。v o r o n o i 图构造研究的问题为如何将给定的平面几何区域划分成一系列有限 个v o r o n o i 几何基元，单个v o r o n o i 几何基元称之为v o r o n o i 网格，因此，v o r o n o i 图的构造也被称为v o r o n o i 网格划分。v o r o n o i 图构造技术在诸多领域都有广泛的 应用，如计算机几何图形学、计算机视觉、科学可视化、g i s 地理信息系统，移 动通信，机器人技术等领域，其研究关注者众多，但由于在各学科领域研究的内 容、层次、角度、需求和应用不尽相同，所采用的方法和面临的挑战也互有差异。 1 1v o r o n o i 图的基本概念 1 1 1v o r o n oi 图的定义 计算几何中定义的v o r o n o i 图的构造是将平面上n 个离散点分为一系列小区 域，每一个离散点代表一个区域，该离散点所在的平面区域是指平面上任意一个 坐标点位置到该离散点的直线距离小于或等于到其它n 1 个离散点直线距离( 等 于时为v o r o n o i 网格边上的点) 。 令平面区域a 上有离散的点集合p = p i , p 2 ，p n ) ，那么v ( v i ) 定义为所有到 点p i 直线距离最小的点的集合： 第1 章绪论 v ( p i ) = p l d i s t a n c e ( p ，p i ) d i s t a n c e ( p ，p j ) ，p a ，j i ，j = 1 ，2 ，”，n 平面区域a 的v o r o n o i 图v ( p ) 为所有p i 点的v o r o n o i 图的集合： v ( p ) = v ( p 1 ) ，v ( p 2 ) ，v ( p n ) ) 为了在限定条件下生成理想的v o r o n o i 网格，在此先引入关于v o r o n o i 网格 的两个基本概念： ( 1 ) v o r o n o i 边：每一个v o r o n o i 网格上的每一条连接线段。 ( 2 ) v o r o n o i 点：三条v o r o n o i 边相交的点( 假设无四点共圆) 。 其中v ( p o 称为一个v o r o n o i 网格。平面区域给定的离散点集如下图1 1 左边 所示，根据离散点集的坐标位置信息划分的平面区域v o r o n o i 网格图见图1 1 右 边所示。 图1 1 平面离散点集v o r o n o i 图的构造 1 1 2v o r o n o i 图主要生成算法 到上个世纪6 0 年代，v o r o n o i 图的概念在自然科学和社会科学领域中已经得 到广泛应用。但因为当时只能依赖入力手工画图，如使用圆规和直尺进行v o r o n o i 图的构造，其缓慢的工作量不能满足广泛的实际应用需求。但随着计算机科学领 域的迅速发展和v o r o n o i 图的大量实际应用需求促使了研究人员成功地设计了多 种关于v o r o n o i 图的构造算法。目前，v o m n o i 生成方法按照构造方式主要分为直 接法和间接法两类：直接法和间接法。直接法即指根据v o r o n o i 的数学定义直接 生成v o r o n o i 图，主要实现算法有半平面的交、增量算法、分治法、平面扫描算 法和离散构造算法，栅格算法等；间接法则利用d e l a u n a y = 角网格是v o r o n o i 网格 的对偶性质，先生成d e l a u n a y = 角剖分，然后构造v o r o n o i 图。以下简要介绍目前 几种常用的直接构造v o r o n o i 图算法的基本设计思路。 ( 1 ) 半平面的交 2 第1 章绪论 利用v o r o n o i 图的定义等式v ( p i ) = n h ( p i ，p j ) 构造1 3 一1 个半平面的交，得到点 p i 的v o r o n o i 图，然后逐点构造其它个点的v o r o n o i 网格得到最终的离散点集p 的 v o r o n o i 图。该算法的时间复杂度n o ( n 2 ) 。 ( 2 ) 增量算法 该算法的基本思想是增量生成的“追加引进”，逐个增加点集中离散点，每 次添加后即时v o r o n o i 图的边。增量算法分为联机增量算法和脱机增量算法，其 中所谓脱机增量算法是指在执行生成算法之前，给出了点集中所有点的位置信 息。而联机增量算法是指，当阶段计算完成时算法才能接受( 请求) 一个新的输入， 或者数据按某种规律到达( 时间间隔) ，这时要求算法在下一个数据到达之前完成 阶段计算。增量算法时间复杂度为：最好情况是o ( n ) ，最坏情况是o ( n 2 ) ，平均 情况是o ( r a o g n ) 。 ( 3 ) 分治法 分治法是目前认为的最优v o r o n o i 图生成算法，其基本思想是：将平面离散 点集p 中的所有的点按x 坐标由小到大排序，然后将点集p 平均划分为p l 和p 2 两部 分，使得p _ p lup 2 。如果点集合p i 和p 2 中含点数目大于4 ，则继续平分点集，直 至子点集的点集。每个子点集各自构造v o m n o i 图后，不断合并相邻子点集，直 至最终得到整个平面离散点集的v ( p ) 。该算法的时间复杂度为o ( n l o g n ) 。 ( 4 ) 平面扫描算法 平面扫描算法通过平面上的一条扫描线由左至右进行扫描。对平面上已扫描 过的部分能够构造出v o r o n o i 图，直到扫描完平面上所有离散点集。该算法的时 间复杂度为o ( r l l o g n ) 。 ( 5 ) 离散构造算法 离散构造算法也称颜色法，给平面点集中每个点分别指定不同的颜色并以其 为圆心，匀速的向外扩张画圆至整个平面，不同颜色的区域相遇产生的交线图即 为v o r o n o i 图。该算法是很直观的一种方法，不仅适用于一般v o r o n o i 图的构造， 也适用于具有障碍物的v o r o n o i 图的构造。 ( 6 ) 栅格算法 该算法利用距离变换，将欧氏空间中两点的距离变换为栅格空间中两点的距 离，进而获取邻元关系，从而得到整个平面区域点集v o r o n o i 图。 3 第1 章绪论 1 2 并行v o r o n oi 网格简介 众所周知，随着近年来计算机科学技术的迅猛发展，高性能计算对计算能力 的需求永无止境，待处理分析的计算数据量越来越大，仅提高单核芯片的运算速 度已经无法满足大规模数据的实时计算处理。同时，由于微处理器工艺、散热等 物理技术的诸多限制，符合摩尔定律发展的处理器已经由原有的提高主频的方式 转变为在芯片内部集成更多处理单元的方式，处理器的发展使得单核处理器的纵 向挖深转变到横向发展的多核处理器时代，多核处理器的概念由此而生。 多核处理器是指在一枚处理器中集成两个或多个完整的内核计算引擎。多核 处理器按照处理器的核结构类型分成同构多核和异构多核两类，同构多核是指在 多核处理器内，集成各个核处理器的结构是相同的，各个核的地位作用平等一致， 而异构多核是指在多核处理器内各个核的结构地位是不同的，根据功能作用分为 主核、副核。因此，随着高性能计算需求的日益增长，多核处理器在高性能计算 中间得到了广泛的普及。根据高性能计算的基本概念：并行程序的加速比决定于 串行部分的并行优化的性能和效率，所以，从理论上分析，异构多核微处理器的 结构显然具有更好的并行加速性能。 1 2 1 并行v o r o n o i 网格算法分类 在目前现有研究中，构造并行v o r o n o i 网格生成算法几乎都是通过区域分解 的方式将问题划分成子问题，各个子问题之间的几何相邻边界称为交界面 ( i m e f f a c e ) 。根据交界面网格生成和子问题内部网格生成之间的先后顺序，并行 v o r o n o i 网格生成算法可分为前置处理类、后置处理类和同步处理类三个类型 【3 2 】： ( 1 ) 前置处理类：先生成子问题之间的交界面网格，后生成子问题内部网格。 ( 2 ) 后置处理类：先生成子问题内部网格，后生成子问题之间的交界面网格。 ( 3 ) 同步处理类：同时生成子问题之间的交界面网格和子问题内部网格。 1 2 2 并行v o r o n o i 网格算法目标 数据并行算法是分而治之( d i v i d ea n dc o n q u e o 思想在并行算法设计上的典型 应用之一，其关键环节包括数据分解和数据映射。数据分解是指将问题分解成子 问题，大规模数据逐步向下划分为小规模数据；数据映射是指将子问题分配到不 4 第1 章绪论 同的处理器上同时计算处理。为提高整个并行算法计算过程的效率，在讨论并行 网格生成算法的效率时，还要关注负载平衡( 1 0 a d i n gb a l a n c e ) 和通信最小化 ( m i n i m i z a t i o no f c o m m u n i c a t i o n s ) ，即分配给每个参与计算的处理器的计算量应近 似相等，且各处理器之间的通信数据应尽可能最小，均衡考虑才能达到并行计算 性能最大化。 目前，已经存在的构造v o r o n o i 图的并行方案有两种方式： ( 1 ) 使用与平面离散点集同等数量级的处理器数目，i i p p = o ( n ) ，复杂度能达 虱j o ( t g n 。l g n t g n ) ，适用于对实时处理要求非常高的场合，但是这样的处理环境 成本也是非常高。 ( 2 ) 使用常数p 个处理器来对任意数量级规模的计算问题进行处理，此类算 法易于实现，算法提速明显，很适用于一般场合的计算。 基于实际情况和算法普适性的考虑，本文使用了第二种并行算法实现策略， 即使用p 个处理器来对v o r o n o i 并行构造问题进行分析。 1 2 3 并行v o r o n o i 网格算法标准 对于衡量一个并行v o r o n o i 网格构造算法的优劣，c h r i s o e h o i d e s 2 5 给出了以 下三个衡量标准： ( 1 ) 稳定l 生( s t a b i l i t y ) ：并行算法是否能获得和串行化算法一样的网格质量。 ( 2 ) 代码复用性( c o d e r e u s i n g ) ：并行算法是否能复用串行化网格生成算法。 ( 3 ) 可扩展性( s c a l a b i l i t y ) - 随着计算规模的扩大，使用同比增长的计算资源， 并行算法是否能保持运行效率的稳定。 本文基于复杂限定情况和不规则多边形区域的情况下生成合理的v o r o n o i 图，离散整个求解区域为一系列的v o r o n o i 网格，之后对模拟问题数学建模，进 行数值分析，将连续方程离散化，而求出每个v o r o n o i 网格质心( 即平面上生成的 离散点集中的点) 的特征量( 如压力，流速等) ，用此特征量的值来代表整个所在 v o r o n o i 网格区域的值，一般在视觉上主要采取用不同的颜色来标识，给人以主 观上的可视直观性。因此，受限v o r o n o i 图的网格质量生成好坏同样取决于在现 实应用中，达到人主观观察上的视觉良好性的结果，而使用到的最少计算量。 s 第1 章绪论 1 3受限v o r o n oi 网格简介 “受限”是指在给定的需划分的平面几何区域中出现了限定情况，如给定点 必须为v o r o n o i 网格节点，v o r o n o i 网格边必须沿给定线段分布，不允许出现跨过 给定线段的v o r o n o i 网格，可见，受限v o r o n o i 网格是指将给定的具有复杂情况的 平面几何区域划分为一系列合理布局的v o r o n o i 网格，从而达到符合实际应用的 需求。 以实际应用中的油藏模拟研究为例，该技术作为油田开发科学计策的工具， 在石油工业中对油、气田的开发方面有着重要的实际应用，而在开发前期的网格 生成技术在油藏数值模拟中占据着重要的地位，通过将油藏模拟区域划分为一系 列合理的v o r o n o i 网格进行区域离散化，将代表性物体性质( 如压力，流速等) 赋值给各个v o r o n o i 网格，来达到模拟物体过程的动态效果。 1 4v o r o n oi 网格在油藏模拟中的应用 在油藏模拟中，应用的网格类型按照网格结构类型可划分为结构网格和非结 构网格，结构网格的定义为在任何网格相交点处的网格数目是一致的，而非结构 在网格相交点处的网格数目是不确定的。最早出现并且应用最为广泛的结构网格 生成技术为笛卡儿坐标系网格，笛卡儿坐标系网格具有生成速度快，数据结构简 单，对后期计算处理简单方便等优点，但其适用范围有限，只能适用形状规则的 平面图形，在描述油臧模拟中复杂地质条件时，不能精确描述边界形状，且存在 严重的网格取向问题。而结构网格一角点网格虽然满足了可以对复杂边界进行描 述，但又失去了网格的正交性。h e i n e m a n 于1 9 8 9 年提出了非结构网格一 p e b i ( p e r p e n d i c u l a rb i s e c t i o n ) 霞j 格来替代在地质模拟领域使用的矩形结构网格。 p e b i 网格指这样的一组网格，任意两个相邻网格块的交界面一定垂直平分相应 网格节点的连线，每个p e b i 节点对应一个p e b i 多边形，p e b i 网格实际上就是 将计算几何中的v o r o n o i 网格概念应用在油藏模拟领域，v o r o n o i 网格加上油藏 模拟领域中可能用到的一些复杂的限定情况( 如不规则边界，断层，井点等) 就 形成了p e b i 网格，即p e b i 网格为最终的根据复杂受限情况生成的受限v o r o n o i 网格。相对于前面所提的结构网格，构造的受限v o r o n o i 非结构网格可以逼近任 意油藏区域形状，同时便于局部加密处理，具有易于描述断层等特殊的限定条件 的优势，受限v o r o n o i 非结构网格目的是尝试建立相对于结构网格更为精细的地 6 第1 章绪论 质模型，提高油藏模拟过程中的计算速度和精度。在油藏模拟中卅现的备个网格 如下图12 所示，其中区域中各种不同的颜色代表了各个网格具有小同的物理性 质。 鳓专鹪 蚓i2 依次为结构网格，非结构网格，角点哪格，受限v o r o n o i ( p e b ) 同格 为了对油减模拟等相关工程领域的某些问题进行更精确的描述和分析，往往 需要大量的计算和生成太规模的受限v o r o n o i 网格数据，而在现有的单机单核的 情况下，生成大规模的网格数据会在计算时间和计算效率上形成严重的性能瓶 颈，并行v o r o n o i 网格生成技术正在在这背景下应运而生。 15目前国内外研究现状和发展趋势 目前，国外的s t a n f o r d 大学在对油减数值模拟技术的研究上处于领先地位， v e r m a 6 提出了二维平面区域下p e b i 网格生成方法，它假定内部边界( 如断层) 和外部边界是直线段。k a p p a 6 考虑了边界附近节点、角点和竖直井周围节点等 限定情况互相之间的干扰，总结了二维的p e b i 网格生成方法。而咀进行地质模 拟研究而著称的s p e 上，只有零星的论文出现在s p e 年会上。 在国内油藏模拟的网格技术研究中，蔡强对采用控制圆算法束实现复杂 条件下生成二维p e b i 网格。安永生3 1 提出了综合运用径向嘲格、笛卡尔网格和 六边形p e b i 网格的混合p e b 网格划分。刘立明1 1 1 提出了综合运用径向网格和 p e b i 刚格的混合网格分别模拟油井区和油减区。谢海兵f 1 2 1 研究了二维的p e b i 网格生成，但是考虑的条件比较简单，只有外部边界和井点，没有考虑其它限定 情况。向祖平 7 】针对油臧任意约束平面多边形区域提出了一种实用的局部j 下交 化网格( p e b d ，士成方法。杨权一 1 卅引入控制圆和桥边解决了维平面复杂断 层限定，对模拟区域进行了更细致的划分。 在p e b i 网格技术应用方面，随着油减数值模拟理论的不断完善和计算机性 第1 章绪论 能的大幅提高，欧美国家开发了相当多的关于地质、地震与油藏模拟类软件，并 逐步趋向成熟。但是由于区域中复杂的限定条件，使得这些软件在p e b i 网格生 成技术上都不够完善。 图1 3 某油藏数值模拟软件生成的地质网格图 图1 3 为某油藏数值模拟软件中所描述生成的地质实例，根据地质数据操作 生成图1 3 中的不规则网格v o r o n o i 图，由图1 3 中可见，在红色水平线左端有 v o r o n o i 网格超出不规则边界。 因此，根据市场现有数值模拟软件和网格生成技术的发展趋势，本文提出一 个成熟的受限v o r o n o i 网格生成软件应具备下列相关性能： 算法应当具有较高的自动性，减少人工干预。 软件在输入数据方面简单方便，与多种通用建模软件具有文件接口。 算法具有较高的通用性和可靠性。 在生成速度上具有高效性。 网格质量要满足计算需要。 数据结构简单便于维护和改进。 8 玲 轮 治 治 玲 珩 船 o 抢 第1 章绪论 1 6论文的研究目的、内容和组织 作为油藏模拟技术的重要一环，受限v o r o n o i 网格生成技术研究在网格生成 质量和网格生成效率面临巨大挑战，也成为了制约提高分析效率的瓶颈之一。因 此，首先要满足生成的符合限定情况和不规则边界的高质量网格，这是基本要求 之一。同时，在现有的单机单核计算能力下，处理大规模数量的网格生成，算法 效率低下，计算时间无法接受，而在大型并行处理机进行数据处理，价格昂贵， 用户操作不便。随着硬件处理器的快速发展，以异构多核为特征的众核处理器已 成为高性能计算的主流趋势，因此，在异构多核环境下构造高质量v o r o n o i 网格 并行化算法的研究有着很重大的意义。 本文的具体组织和安排如下： 第一章讨论了本课题的研究v o r o n o i 网格的应用背景及相关基本定义概念， 阐述了v o r o n o i 网格并行化算法分类、目标，以及受限v o r o n o i 网格在油藏模拟 实例的应用优势等相关内容。 第二章阐述了利用d e l a u n a y 网格与v o r o n o i 网格之间的直线对偶性质，基于 间接法先生成d e l a u n a y 网格，然后连接各个相邻d e l a u n a y 网格的外接圆心来转 换成相应的v o r o n o i 网格，并给出实验结果，由于算法最后对d e l a u n a y 网格进行 了优化处理，生成的网格模型质量较高，但运行速度较慢。 第三章重点阐述了基于直接法分治算法生成受限v o r o n o i 网格，详细分析了 如何在各种限定情况的干扰下合理布点，并实现了并行布点划分方法，利用分治 法是最优算法且易于并行化的优势，构造双向循环链表对分治法进行了并行化， 最后对边界外的点做了特殊处理，使用逐步并行的并行方案来是实现整个算法并 行一体化。 第四章总结了本文提出的限v o r o n o i 图并行生成算法，通过在大型惠普并行 机和异构多核机器上对于分治法并行化算法做了实验结果分析，最后将实现的算 法与当前的一些受限v o r o n o i 图算法在计算时间和网格生成质量上进行了分析比 较。 第五章对基于v o r o n o i 图技术未来的研究方向，进行了展望。 基于本文的研究内容，可以把该文的研究成果应用于模拟流体流动、建模的 绘制、各种数值模拟等领域，同时，结合异构多核处理器和并行计算高速发展的 时代特点，用于各种模拟问题的实时计算。 9 第2 章基于间接法实现v o r o n o i 图 第二章基于间接法实现v o r o n o i 网格 目前，构造v o r o n o i 图按照实现方式有两种，直接实现v o r o n o i 图和间接实 现v o r o n o i 图。直接构造v o r o n o i 图有若干种方法，已在第一章做了引述。本章 根据d e l a u n a y 图是v o r o n o i 图的直线对偶图特性，先实现d e l a u n a y 图的构造， 再连接已构造好的相邻d e l a u n a y 三角形网格的外接圆心转换成相应的v o r o n o i 图，如下图2 1 所示。 一 d e l 蛳a y 图连线 一一 v o r o x w i 图连线 l ， 图2 1d e l a u n a y 图与v o r o n o i 图的转换 2 1d e i a u n a y 三角剖分的基本概念 2 1 1d e i a u n a y 三角剖分定义 对平面区域上离散的点集合p = p l ，p 2 ，p 。 ，令生成的d e l a u n a y 三角网格边 集合为t r i ( p ) ，若t r i ( p ) 满足条件：设三角网格的边e i j = ( p i ，p j ) ，对于三角网格 p i p i p k ，e i j ，e i k ，e k j 属于t r i ( p ) ，则该三角网格p i p j p k 的外接圆中不能包含平面 离散点集合p 中的任何其它离散点，即四点不能共圆原则。那么，则称t r i ( p ) 是 平面区域离散点集p 的d e l a u n a y 三角剖分。根据d e l a u n a y 图是v o r o n o i 图的直 线对偶图特性，依次连接t r i ( p ) 中所有相邻三角形的外接圆心，则转换为相应的 v o r o n o i 图。如下图2 2 ，由左及右，平面上的离散点集，三角剖分实现，再转 换v o r o n o i 对偶网格，显示了利用三角剖分间接实现v o r o n o i 网格的划分构成。 1 0 第2 章基于间接法实现v o r o n o i 图 图2 2d e l a u n a y 三角剖分定义 d e l a u n a y 三角剖分必须符合两个重要的准n 2 4 ，即空圆准则和最大化最小 角准则： ( 1 ) 空圆准则：平面区域离散点集的三角网格划分中任意四点不能共圆，即 在d e l a u n a y 三角形剖分图中任一个三角形的外接圆区域内不会存在平面离散点 集中的其它点。 ( 2 ) 最大化最小角准则：在d e l a u n a y 三角形剖分图中任一个形成的三角形最 小角最大，是指在相互交换任意两个相邻的三角形构成凸四边形的对角线后，六 个内角的最小角不再变大。 2 1 2d e l a u n a y 三角剖分性质 d e l a u n a y 三角剖分具有最接近、唯一性、最优性、最规则、区域性等特点， 具体数学描述如下： ( 1 ) 生成d e l a u n a y 图上所有离散点都在一个平面上，生成的d e l a u n a y 图是 一个平面图。 ( 2 ) 由四点共圆特性将d e l a u n a y 三角剖分化为两类，有四点共圆平面离散点