CN102662641A - 一种基于cuda的种子分布数据并行获取方法 - Google Patents

Abstract

一种基于CUDA的种子分布数据并行获取方法，该方法在CUDA架构下，在CPU上执行的主机端和在GPU上执行的设备端；所述种子分布数据并行获取方法包括以下步骤：(1)、输入数据内存空间分配并初始化；(2)、输入数据显存空间分配；(3)、内存到显存的输入数据传输；(4)、输出数据显存空间分配；(5)、输出数据内存空间分配；(6)、设置设备端内核执行配置，调用设备端内核函数计算种子分布。本发明提供一种在确保精度的同时获得较好的计算性能、降低演化仿真时间的基于CUDA的种子分布数据并行获取方法。

Description

一种基于CUDA的种子分布数据并行获取方法

技术领域

本发明涉及一种种子分布数据并行获取方法。

背景技术

种子分布计算是森林群落更新繁殖体的主要来源，是森林演化计算的重要子过程。森林中种子的分布情况是成熟成年树生产种子并扩散的结果。种子通过扩散进入森林样地，随着土壤基质的促进作用萌发成幼苗，幼苗受环境条件影响生长成幼年树甚至成年树。生长状况不好的树木会出现衰弱死亡，而个体成熟的成年树继续生产种子，为森林提供天然更新的物质基础。种子分布、萌芽、树木生长和死亡构成森林演化的动态循环过程，反映森林群落的发展情况。

种子分布主要受个体胸径(高1.3米处的树干直径)大小和树种两个因素限制，胸径越大表示越成熟，种子生产能力越强；种子的扩散与树种特性有关，表示为扩散内径曲线。一般来说，某种子分布点距离种子源(产生种子的成年树)越近，扩散到该点的种子数越多；距离越远，数量越小。不同树种扩散距离不同。森林样地里某个种子分布点的种子总量受所有树木个体生产扩散累积。

通常，森林样地被简化模拟成矩形区域，由X-Y坐标系统离散化为多个面积相同的样地单元，以每个样地单元的中心点表示样地单元的种子分布点。本系统采用的种子分布模型如下所示，示意图见图1：

$R_i=\frac{1}{\mathrm{\η}}\mathrm{STR}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\left({\left(\frac{{\mathrm{DBH}}_k}{30}\right)}^{\mathrm{\β}}{e}^{-{\mathrm{Dd}}_{\mathrm{ik}}^{\mathrm{\θ}}}\right)$

其中，R_i表示第i块样地单元内分布的种子数，N表示成年树数量，DBH_k表示第k棵成年树的胸径，d_ik表示第k棵成年树到样地单元i中心点的距离。其余字符均表示具有生物意义的常量参数。因此完成所有样地单元的种子分布计算包含两层循环：遍历所有成年树木个体，计算这些个体分布到某个样地单元的种子总数；重复以上步骤，完成每个样地单元的计算。

然而森林是个复杂庞大的系统，不但样地分布广，而且数量繁多，精确计算样地内各样地单元的种子分布数需要耗费大量时间。根据实验结果得到在1000m×500m的森林样地上，若每平方米代表一个样地单元，模拟100000棵树的种子分布计算花费146小时。

目前对森林演化计算领域的加速优化方法研究较少。杜克大学Govindarajan等人提出种子分布的聚簇优化算法。该算法依赖于种子扩散距离关系。根据种子扩散内径曲线分析，若某成年树距离种子分布点较远，其空间位置的微小变化对此成年树扩散到该点的种子数影响不大。利用此空间距离关系，通过预处理方式将与某分布点超过一定距离且空间分布密集的成年树群聚合为超级树，以超级树对该收集点的种子贡献数取代成年树群中个体的分别计算，减少总的计算次数。该算法以牺牲精度获取计算效率的提高。然在样地规模较大时，预处理部分较耗时，且当种子扩散距离依赖关系不明显时不再适用。

发明内容

为了克服已有种子分布数据获取技术的不能兼顾获取精度和计算性能、耗时较大的不足，本发明提供一种在确保获取精度的同时获得较好的计算性能、降低演化仿真时间的基于CUDA的种子分布数据并行获取方法。

本发明解决其技术问题所采用的技术方案是：

一种基于CUDA的种子分布数据并行获取方法，该方法在CUDA架构下，在CPU上执行的主机端和在GPU上执行的设备端；所述种子分布数据并行获取方法包括以下步骤：

(1)、输入数据内存空间分配并初始化：根据样地单元规模分配内存空间存放每个样地单元中心坐标；根据成年树规模分别分配内存空间存放成年树的胸径、空间坐标并初始化；根据种子分布模型所需常量参数分配内存空间并初始化；

(2)、输入数据显存空间分配：调用CUDAAPI函数cudaMalloc在显存的全局存储器中分配空间用于存放样地单元和成年树数据；

(3)、内存到显存的输入数据传输：调用CUDA API函数cudaMemcpy，将以上输入数据从内存复制到全局存储器对应区域，并调用cudaMemcpyToSymbol将各常量参数直接复制到常量存储器；

(4)、输出数据显存空间分配：根据内核配置情况分配显存空间，用于存放设备端的计算结果；

(5)、输出数据内存空间分配；

(6)、设置设备端内核执行配置，调用设备端内核函数计算种子分布，设备端内核函数计算步骤如下：

(6.1)根据CUDA内建函数获得线程块的Y索引，即目标样地单元，记为bid_y；

(6.2)在0号线程中根据bid_y读取目标样地单元的空间坐标保存到共享存储器；初始化数组seeds_thread各元素为0；调用CUDA同步函数__syncthreads()，使得块内所有线程在初始化完成后再进行种子分布数计算累加；

(6.3)根据CUDA内建函数获得该线程所在线程网格的X索引，记为tid_x，表示参与该线程计算的成年树集；设置寄存器变量，对每棵成年树i，i＜n(n值为常量)，根据索引tid_x×n+i获得成年树的胸径和空间位置，按照种子分布计算公式计算种子数并累加到寄存器变量；待n次迭加完成后，将该寄存器变量写入共享存储器seeds_thread数组中；

(6.4)调用CUDA同步函数__syncthreads()，采用树形规约算法，通过多线程并发执行实现每个线程的种子数累加操作，得到该线程块的种子贡献总数；

树形规约算法中，s表示数组长度，初始值为T，对于seeds_thread数组的前半部分元素，每一个与比其索引大s的元素相加后更新其值，每次迭代后将s值减半；

(6.5)在0号线程根据线程块索引将结果保存至全局存储器seeds_block中，得到种子分布。

进一步，所述种子分布数据并行获取方法还包括以下步骤：

(7)、显存到内存的输出数据传输，调用CUDA API函数cudaMemcpy，将设备端的种子分布计算结果从显存复制到内存相应区域；

(8)处理计算结果；

(9)释放显存空间，释放内存空间。

再进一步，所述步骤(6)中的内核执行配置包括根据种子分布的计算特性划分设备端线程网格和存储器访问模式设置；

线程网格用二维表示，X维度表示贡献种子的成年树个体，Y维度表示分布到种子的各样地单元，每个线程块计算一定数量成年树对一个样地单元的种子贡献数；相同Y维度的线程块负责对同一个样地单元的种子数计算；同时；在全局存储器中分配一维线性数组seeds_block用于存储每个线程块贡献的种子数；该数组中属于同一样地单元的元素值之和即为该样地单元的种子数终值；

假设共有A棵成年树，每个线程块包含T个线程，每个线程处理n棵成年树的计算，那么X维度共划分为A/(n×T)个线程块；Y维度的线程块数取决于样地单元总数和显卡配置；若以TL表示该内核执行的总线程数上限，那么Y维度的线程块数为(TL×n)/A；若样地单元总数大于该数值，那么样地单元被分成若干批执行，通过多次执行内核函数完成对所有样地单元的计算；

在所述内核函数中，利用共享存储器存储两部分数据，一是目标样地单元数据；二是每个线程n棵成年树贡献的种子数和，表示为float型数组seeds_thread，该内核共需开辟的共享存储器大小为(T+2)×4字节，其中T表示seeds_thread数组长度，即块内线程数；2表示目标样地单元X，Y方向上的空间坐标；4表示每个float型占4个字节。

更进一步，所述步骤(6)的内核函数执行完毕后，按照步骤(7)将seeds_block数组拷贝回内存。根据所述步骤(8)，在CPU端将seeds_block中属于同一样地单元的块结果再做累加操作，获得各样地单元的种子数终值。

本发明的技术构思为：通过分析发现，种子分布计算公式逻辑简单，计算挑战主要在于参与计算的大规模森林场景数据。根据其计算密集度高、数据可高度并行化的计算特性，可以将总的计算任务通过一定方式分解为细粒度的子任务并行处理。由于两个种子分布点间的种子数计算具有独立性，且不同成年树产生的种子分布到某一点的过程也没有关联性，种子分布可被分解为多个独立、计算任务一致且可并行的计算单元。每一计算单元处理一定数量成年树分布到某一收集点的种子数。

传统的CPU处理器体系结构对数据并行应用支持不佳，CPU上更多的晶体管是用于控制分支和cache，而非实际的运算工作，无法提供强大的计算能力。与CPU相反，GPU拥有丰富的计算资源和简单的控制逻辑资源，在计算能力和存储器带宽上有明显的优势，具有支持数据并行应用的基础。算术密度越高，数据并行度越好的计算任务越适合GPU处理。近年来，随着图形硬件的不断发展和可编程性的提高，GPU除了被用于图形处理以外，也吸引了越来越多研究者将重点转移到解决通用计算任务上来。尤其是NVIDIA公司计算统一设备架构(CUDA)的发布，则让GPU在通用计算领域真正大量应用，其应用广泛分布于物理模拟、科学计算、信号处理、媒体处理等领域。

本发明的有益效果主要表现在：(1)思路新颖，非经典方法。突破传统基于数据结构的CPU加速技术，具有较大创新意义。

(2)开发环境配置简单，开发难度低。CUDA采用了统一处理架构，将GPU上的计算功能和图形处理功能完全区分开，可以更加有效地利用GPU进行通用计算。CUDA同时也提供了让开发者能够使用C作为高级编程语言的软件环境，为编程人员提供从CPU到GPU的平稳过渡模式而不必重新学习语法。

(3)性价比高。一块高端的NVIDI显卡售价几千元，CUDA能够利用NVIDIA的GPU并行计算引擎高效地解决复杂的计算任务。相较之下，利用集群技术在获得相同加速效率的情况下花费的成本较高，性价比较CUDA低。

(4)计算精度高。可靠的计算精度一直是GPU设计者不断追求的目标，通过增加对浮点数据的处理功能，提高渲染精度。支持浮点运算让GPU具备通用计算的基础，基本浮点格式对应C/C++类型中的float和double。对于本应用中种子分布计算，单精度浮点型已满足精度需求。实验结果显示CPU和GPU计算结果的误差率仅占10-8，计算精度较高。

(5)适用性强。CUDA并行体系非常适合处计算密集度高、数据并行化好的问题。种子分布的特点是复杂流控的要求非常少，同一个模型公式在多个数据上执行，符合CUDA并行体系的计算要求。目前生态学家对林木本身的生物学特性和外界环境对其影响过程的研究不断深入，种子分布模型公式也在不断细化，但基本结构并未改变，仍支持同一模型公式在多数据上并行执行的特性。以本申请设计的CUDA线程网格结构为参考框架，根据实际的生长方程适当调整内核计算方法，可适用于不同的种子分布模型。

附图说明

图1是种子分布模型示意图。

图2是种子分布线程网格划分适宜图。

图3是存储器访问设计示意图。

图4是树形规约算法示意图。

具体实施方式

下面结合附图对本发明作进一步描述。

参照图1～图4，一种基于CUDA的种子分布数据并行获取方法，在CUDA架构下，该方法分为两个部分：在CPU上执行的主机端和在GPU上执行的设备端；主机端程序包括以下步骤：

(1)、输入数据内存空间分配并初始化。根据样地单元规模分配内存空间存放每个样地单元中心坐标；根据成年树规模分别分配内存空间存放成年树的胸径、空间坐标并初始化；根据种子分布模型所需常量参数分配内存空间并初始化。

(2)、输入数据显存空间分配。调用CUDAAPI函数cudaMalloc在显存的全局存储器中分配空间用于存放样地单元和成年树数据。

(3)、内存到显存的输入数据传输。调用CUDA API函数cudaMemcpy(参数cudaMemcpyKind设置为cudaMemcpyHostToDevice)，将以上输入数据从内存复制到全局存储器对应区域，并调用cudaMemcpyToSymbol将各常量参数直接复制到常量存储器。

(4)、输出数据显存空间分配。根据内核配置情况分配显存空间，用于存放设备端的计算结果。

(5)、输出数据内存空间分配。内存大小与步骤4中显存大小相同。

(6)、设置设备端内核执行配置，调用设备端内核函数计算种子分布，设备端内核函数计算步骤如下：

(6.1)根据CUDA内建函数获得线程块的Y索引，即目标样地单元，记为bid_y；

(6.2)在0号线程中根据bid_y读取目标样地单元的空间坐标保存到共享存储器。初始化数组seeds_thread各元素为0。调用CUDA同步函数__syncthreads()，保证块内所有线程在初始化完成后再进行种子分布数计算累加，避免发生访问出错的情况。

(6.3)根据CUDA内建函数获得该线程所在线程网格的X索引，记为tid_x，表示参与该线程计算的成年树集。设置寄存器变量，对每棵成年树i(i＜n，n值为常量)，根据索引tid_x×n+i获得成年树的胸径和空间位置，按照种子分布计算公式计算种子数并累加到寄存器变量。待n次迭加完成后，将该寄存器变量写入共享存储器seeds_thread数组中。

(6.4)调用CUDA同步函数__syncthreads()，保证完成块内所有线程的上述计算与存储。采用树形规约算法，通过多线程并发执行实现每个线程的种子数累加操作，得到该线程块的种子贡献总数。

树形规约算法实现如图4所示，s表示数组长度，初始值为T。对于seeds_thread数组的前半部分元素，每一个与比其索引大s的元素相加后更新其值，每次迭代后将s值减半。

(6.5)在0号线程根据线程块索引将结果保存至全局存储器seeds_block中，到此完成种子分布的内核计算。

(7)、显存到内存的输出数据传输。调用CUDA API函数cudaMemcpy(参数cudaMemcpyKind设置为cudaMemcpyDeviceToHost)，将设备端的种子分布计算结果从显存复制到内存相应区域。

(8)处理计算结果。

(9)释放显存空间，释放内存空间。

所述步骤(6)中的内核执行配置，主要包括根据种子分布的计算特性划分设备端线程网格和存储器访问模式设置。

线程网格用二维表示，如图2。X维度表示贡献种子的成年树个体，Y维度表示分布到种子的各样地单元。每个线程块计算一定数量成年树对一个样地单元的种子贡献数。相同Y维度的线程块负责对同一个样地单元的种子数计算。同时我们在全局存储器中分配一维线性数组seeds_block用于存储每个线程块贡献的种子数。该数组中属于同一样地单元的元素值之和即为该样地单元的种子数终值。为进一步增加每个线程的计算量以隐藏访问延迟，我们让每个线程处理(n＞＝1)棵成年树的种子贡献数。

假设共有A棵成年树，每个线程块包含T个线程，每个线程处理n棵成年树，那么X维度共划分为A/(n×T)个线程块(若不能整除，则对结果向上取整以保证所有成年树都参与计算)。Y维度的线程块数取决于样地单元总数和显卡配置。若以TL表示该内核执行的总线程数上限(该数据可根据实验获得)，那么Y维度的线程块数为(TL×n)/A。若样地单元总数大于该数值，那么样地单元被分成若干批执行，通过多次执行内核函数完成对所有样地单元的计算。

存储器访问模式设置表示内核函数执行时使用的共享数据大小。在该内核函数中，我们利用共享存储器存储两部分数据，一是目标样地单元数据，如空间坐标，供块内所有线程共享；二是每个线程n棵成年树贡献的种子数和，表示为float型数组seeds_thread。因此，该内核共需开辟的共享存储器大小为(T+2)×4字节，其中T表示seeds_thread数组长度，即块内线程数；2表示目标样地单元X，Y方向上的空间坐标；4表示每个float型占4个字节。

所述步骤(6)的内核函数执行完毕后，按照步骤(7)将seeds_block数组拷贝回内存。根据所述步骤(8)，在CPU端将seeds_block中属于同一样地单元的块结果再做累加操作，获得各样地单元的种子数终值。

Claims (4)

1.一种基于CUDA的种子分布数据并行获取方法，其特征在于：该获取方法在CUDA架构下，在CPU上执行的主机端和在GPU上执行的设备端；所述种子分布数据并行获取方法包括以下步骤：

(1)、输入数据内存空间分配并初始化：根据样地单元规模分配内存空间存放每个样地单元中心坐标；根据成年树规模分别分配内存空间存放成年树的胸径、空间坐标并初始化；根据种子分布模型所需常量参数分配内存空间并初始化；

(2)、输入数据显存空间分配：调用CUDAAPI函数cudaMalloc在显存的全局存储器中分配空间用于存放样地单元和成年树数据；

(3)、内存到显存的输入数据传输：调用CUDA API函数cudaMemcpy，将以上输入数据从内存复制到全局存储器对应区域，并调用cudaMemcpyToSymbol将各常量参数直接复制到常量存储器；

(4)、输出数据显存空间分配：根据内核配置情况分配显存空间，用于存放设备端的计算结果；

(5)、输出数据内存空间分配；

(6)、设置设备端内核执行配置，调用设备端内核函数计算种子分布，设备端内核函数计算步骤如下：

(6.1)根据CUDA内建函数获得线程块的Y索引，即目标样地单元，记为bid_y；

(6.2)在0号线程中根据bid_y读取目标样地单元的空间坐标保存到共享存储器；初始化数组seeds_thread各元素为0；调用CUDA同步函数__syncthreads()，使得块内所有线程在初始化完成后再进行种子分布数计算累加；

(6.3)根据CUDA内建函数获得该线程所在线程网格的X索引，记为tid_x，表示参与该线程计算的成年树集；设置寄存器变量，对每棵成年树i，i＜n(n为常量)，根据索引tid_x×n+i获得成年树的胸径和空间位置，按照种子分布计算公式计算种子数并累加到寄存器变量；待n次迭加完成后，将该寄存器变量写入共享存储器seeds_thread数组中；

(6.4)调用CUDA同步函数__syncthreads()，采用树形规约算法，通过多线程并发执行实现每个线程的种子数累加操作，得到该线程块的种子贡献总数；

树形规约算法中，s表示数组长度，初始值为T，对于seeds_thread数组的前半部分元素，每一个与比其索引大s的元素相加后更新其值，每次迭代后将s值减半；

(6.5)在0号线程根据线程块索引将结果保存至全局存储器seeds_block中，得到种子分布。

2.如权利要求1所述的一种基于CUDA的种子分布数据并行获取方法，其特征在于：所述种子分布数据并行获取方法还包括包括以下步骤：

(7)、显存到内存的输出数据传输，调用CUDA API函数cudaMemcpy，将设备端的种子分布计算结果从显存复制到内存相应区域；

(8)处理计算结果；

(9)释放显存空间，释放内存空间。

3.如权利要求1或2所述的一种基于CUDA的种子分布数据并行获取方法，其特征在于：所述步骤(6)中的内核执行配置包括根据种子分布的计算特性划分设备端线程网格和存储器访问模式设置；

线程网格用二维表示，X维度表示贡献种子的成年树个体，Y维度表示分布到种子的各样地单元，每个线程块计算一定数量成年树对一个样地单元的种子贡献数；相同Y维度的线程块负责对同一个样地单元的种子数计算；同时；在全局存储器中分配一维线性数组seeds_block用于存储每个线程块贡献的种子数；该数组中属于同一样地单元的元素值之和即为该样地单元的种子数终值；

假设共有A棵成年树，每个线程块包含T个线程，每个线程处理n棵成年树，那么X维度共划分为A/(n×T)个线程块；Y维度的线程块数取决于样地单元总数和显卡配置；若以TL表示该内核执行的总线程数上限，那么Y维度的线程块数为(TL×n)/A；若样地单元总数大于该数值，那么样地单元被分成若干批执行，通过多次执行内核函数完成对所有样地单元的计算；

在所述内核函数中，利用共享存储器存储两部分数据，一是目标样地单元数据；二是每个线程n棵成年树贡献的种子数和，表示为float型数组seeds_thread，该内核共需开辟的共享存储器大小为(T+2)×4字节，其中T表示seeds_thread数组长度，即块内线程数；2表示目标样地单元X，Y方向上的空间坐标；4表示每个float型占4个字节。

4.如权利要求2所述的一种基于CUDA的种子分布数据并行获取方法，其特征在于：所述步骤(6)的内核函数执行完毕后，按照步骤(7)将seeds_block数组拷贝回内存。根据所述步骤(8)，在CPU端将seeds_block中属于同一样地单元的块结果再做累加操作，获得各样地单元的种子数终值。

Images