一种大整数乘法Comba算法基于OpenMP的并行实现方法

摘要

本发明公开了一种大整数乘法Comba算法基于OpenMP的并行实现方法，基于64位无符号长整型整数操作，通过添加三个临时数组存储加乘操作计算得到的中间结果，从而解决加乘运算与进位运算的数据相关性，将加乘操作与进位操作分开执行。在加乘操作阶段，基于中间结果每个数位求取时的计算独立性，通过OpenMP多线程编程采用动态调度策略实现加乘操作阶段的并行化，而进位阶段仍然串行执行来并行化Comba算法，提高算法效率。

说明书

技术领域

本发明涉及一种大整数乘法Comba算法基于OpenMP)Open Multi-Processing，一种开 源多核并行技术)的并行实现方法，该方法通过开启多个线程并行计算，提高计算机多核利 用率，从而提高算法效率。

背景技术

大整数乘法是大整数运算的核心运算之一，大多数公钥加密系统算法包括像RSA、ECC、 Diffie-Hellman等中，都要用到大整数乘法，Comba算法则是当大整数数位较小时调用的重 要底层算法。其算法复杂性在于大量的加乘运算以及进位运算，当数据达到一定规模时，将 消耗大量的内存空间和计算时间。Comba算法作为重要的大整数底层乘法运算如果能够提高 运行速度，整个计算的运行效率将会大大改善，在实际应用中有着十分重要的作用。

伴随着多核处理器和加速部件的普及，通过充分挖掘算法本身的并行度将算法并行化， 高效地提升算法性能，已经成为一种行之有效的方法。目前Comba算法优化方面，Vladislav  Kovtun和Andrew Okhrimenko曾通过采用三个64位整型数据存储中间结果以实现延迟进位 (Vladislav Kovtun，Andrew Okhrimenko.Approaches for the performance increasing of software  implementation of integer multiplication in prime fields.[S.l.].2009.)从而改进Comba算法，其 中该实现算法中基为2³²；并基于该改进的算法实现采用OpenMP中section任务分担策略双 线程并行，与OpenMP中for任务分担及静态数据分配策略进行多线程并行两种并行方案) Vladislav Kovtun，Andrew Okhrimenko.Approaches for the parallelization of Software  Implementation of Integer Multipication.[S.l.].2012.)，与单线程程序相比，可分别达到1.5倍及 2倍的加速。其实现主要基于32位单精度整数运算，并行化方法主要基于采用延迟进位技巧 的串行Comba算法，将乘积的计算分成前后两部分数位计算，其section任务分担策略双线 程并行依据乘积的前半部分数位以及后半部分数位计算的独立性，将乘积数位的计算相应地 分成两部分分配给两个线程并行计算，此种方式可实现双线程并行；而其采用for任务分担 策略实现的基于计算资源的多线程并行，首先开启相应数目的线程，采用静态数据分配策略 将前半部分数位计算分配给线程并行计算，计算完任务的线程接着计算乘积的后半部分数 位，但静态数据分配策略容易引起负载不均衡，造成资源的浪费，影响算法性能的提升。

Comba算法基于传统手算乘法，传统的手算乘法在进行运算时将乘数与被乘数的每一位 进行相乘，每相乘一次，便实施进位操作，而Comba算法可首先通过多次加乘操作得到中 间结果的一个数位，即中间结果的每个数位是乘数与被乘数的多对对应位相乘后相加的结 果，获取中间结果的数位后，再向高位进位，然后计算下个数位，从而可顺序求得最终结果。 采用大整数运算领域著名开源库GMP中的大整数数据结构，其定义如下：

其中_mpz_struct即为定义的大整数数据结构，_mp_alloc为大整数动态数组的长度， _mp_size为大整数动态数组已用的长度和符号，_mp_d为指向大整数动态数组首地址的指 针，且大整数每个数位的数据类型为64位无符号长整型整数。

基于上述大整数数据结构，设存储乘数u、被乘数v数组分别为up、vp指向，且乘数与 被乘数位数分别为un和vn位，计算乘积r时先根据如下公式(1)计算其中间结果的第ir数位 r_ir：

$r_{\mathrm{ir}}=\underset{\mathrm{iu}+\mathrm{iv}=\mathrm{ir}}{\Sigma }u{p}_{\mathrm{iu}}*v{p}_{\mathrm{iv}}---\left(1\right)$

算法根据公式(1)循环获得中间结果的每一个数位r_ir(0≤ir<(un+vn))，在求取 r_ir(0≤ir<(un+vn))时，其中乘数下标iv和被乘数下标iu的初始值可按照如下公式设置。

iv←MIN(ir,vn-1)                  (2)

iu←ir-iv                 (3)

Vladislav Kovtun和Andrew Okhrimenko所提供的并行化方法虽然达到了很好的效果， 但是由于其将数位的计算或划分为两个section进行并行，或分为前后两个部分进行for任务 分担并行，且均采用静态调度策略，使得算法并行度低，容易造成负载不均衡而影响算法效 率，而此缺点和不足正是下面技术需要克服的，本发明将提供一种并发度高且负载更为均衡 的并行算法。

发明内容

本发明的解决的问题是：克服现有技术并发度低和负载不均衡的问题，提供一种大整数 乘法Comba算法基于OpenMP的并行实现方法，通过OpenMP多线程编程采用动态调度策 略实现加乘操作阶段的并行化，而进位阶段仍然串行执行来并行化Comba算法，提高算法 效率。

Comba算法的核心操作为加乘运算与进位运算，而进位运算有较强的顺序性，加乘操作 每个数位的计算则具有较强的独立性。为了充分挖掘Comba算法的并行性，本发明基于64 位单精度整数操作，通过解决加乘运算与进位运算的数据相关性，将加乘操作与进位操作分 开，并在加乘操作阶段进行并行优化，进位操作阶段串行执行。

算法的核心操作为64位无符号长整型数的乘加运算，当两个64位无符号整数相乘时， 结果为127或者128位，而目前大部分计算机平台整型数据长度最大为64位，因此需要两 个64位无符号长整型变量分别存储乘积的低64位和高64位。乘积r中间结果的每一个数位 是一对或多对64位无符号长整型数相乘后相加的结果，最终结果会大于或等于127位，所 以中间结果需采用3个64位无符号长整型临时变量存储。如果采用高级语言实现，每进行 一次加乘操作，便需要进行两次进位判断，因此采用内嵌汇编语言的方法实现此核心操作以 消除每次求中间结果的数位时的多次进位判断操作；为了在加乘操作阶段进行并行化，需要 解决加乘操作与进位操作阶段的数据相关性问题，添加三个unsigned long int类型的临时数 组用于存储通过加乘操作求得的中间结果，而在求中间结果的过程中不进行进位操作，待中 间结果求取完毕后再统一进行进位操作；在加乘操作阶段，基于中间结果每个数位求取时的 计算独立性，可通过OpenMP多线程编程采用动态调度策略实现加乘操作阶段的并行化，而 进位阶段仍然串行执行。

本发明的并行方案包括以下方面：

1)参与运算的大整数基为2⁶⁴，申请三个长度为rn的64位无符号长整型数组，分别由 指针lpl，hpl，cl指向，即数组lpl[rn]，cl[rn]以及hpl[rn]，其中rn＝un+vn，un和vn可 能不等，un和vn分别为乘数和被乘数的位数；

2)循环计算乘积的中间结果，缓存在步骤1)申请的三个临时数组中，其中中间结果的 每一位根据公式$\left(\mathrm{lpl}\right[\mathrm{ir}],\mathrm{cl}[\mathrm{ir}],\mathrm{hpl}[\mathrm{ir}\left]\right)=\underset{\mathrm{iu}+\mathrm{iv}=\mathrm{ir}}{\Sigma }u{p}_{\mathrm{iu}}*v{p}_{\mathrm{iv}}(0\le \mathrm{ir}<\mathrm{rn})$计算，lpl[ir]，cl[ir]和hpl[ir] 分别存储中间结果ir位的低64位，中间64位以及高64位，且依据如下方式循 环计算：(rax,rdx)＝up_iu*vp_iv，(lpl[ir],c1)＝rax+lpl[ir]，(cl[ir],c2)＝cl[ir]+hpl[ir]+c1， hpl[ir]＝hpl[ir]+c2，此循环体采用内嵌汇编的方式实现，rax和rdx分别存储up_iu*vp_iv的低 64位和高64位，其中up_iu和vp_iv分别表示乘数u的第iu位和被乘数v第iv位，c1和c2表示进 位；

3)由步骤2)所述$\left(\mathrm{lpl}\right[\mathrm{ir}],\mathrm{cl}[\mathrm{ir}],\mathrm{hpl}[\mathrm{ir}\left]\right)=\underset{\mathrm{iu}+\mathrm{iv}=\mathrm{ir}}{\Sigma }u{p}_{\mathrm{iu}}*v{p}_{\mathrm{iv}}(0\le \mathrm{ir}<\mathrm{rn})$的循环计算过程采用 OpenMP的for任务分担策略进行并行化，开启多个计算线程并行计算，每个线程每次计算 由步骤2)所述的中间结果一个数位的计算流程，线程数可根据当前计算平台可用计算资源 决定，且给每个线程分配计算任务即中间结果数位的计算时采用动态调度策略，其中任务分 割粒度根据开启的线程数目以及计算机资源通过实验获得。本实验中将中间结果的数位分配 给8线程并行求取，数据分割粒度为rn/64，即将数位计算任务共分为组，8 个线程并行执行，每个线程每次领取rn/64个数位计算任务，率先执行完任务的线程接着领 取下组任务，无需等待；

4)对步骤2)和步骤3)所得的lpl[rn]，cl[rn]以及hpl[rn]三个数组中存储的中间结果 通过进位操作将最终结果存入rp[rn]数组中，设rlpl，rcl和rhpl存储进位数，则 rp[ir](0≤ir<rn)以及进位依据如下方式计算：(rlpl,c1)＝rlpl+lpl[ir]， (rcl,c2)＝rcl+cl[ir]+c1，rhpl＝rhpl+hpl[ir]+c2，rp[ir]＝rlpl，rlpl＝rcl，rcl＝rhpl， rhpl＝0，其中(rlpl,c1)＝rlpl+lpl[ir]，(rcl,c2)＝rcl+cl[ir]+c1，rhpl＝rhpl+hpl[ir]+c2， 此三个计算步骤采用内嵌汇编的方式实现，c1和c2表示进位。

本发明的有益效果：

本发明基于64位单精度整数运算，通过添加三个临时数组存储乘积的中间结果，将加 乘操作与进位操作分开执行，从而降低了数据依赖性，通过多线程并行计算，显著提高算法 执行效率。

与上述Vladislav Kovtun和Andrew Okhrimenko实现的并行化方法相比，本发明通过添 加临时数组存储中间结果的方式实现了每个数位计算的独立性，数位的计算并不需要分成前 后两个阶段进行，进而可以实现基于计算资源的多线程并行；并行化时采用动态数据分配策 略，每次分配给线程一定粒度的数据计算任务，其中粒度采用实验的方式确定，更有利于负 载均衡，提升效率，可以达到更优的效果。

附图说明

图1是Comba算法并行策略的程序流程图；

图2是取rn＝36时采用静态数据分配策略的线程任务分配图；

图3是取rn＝36时采用动态数据分配策略一种可能的线程任务分配图；

图4是串行Comba算法与并行Comba算法运行结果对比柱状图。

具体实施方式

下面结合附图及实施例对本发明进行详细说明。

如图1所示，本发明具体实现如下：

1)输入基为2⁶⁴乘数u和被乘数v，其位数分别为un和vn，un和vn可不等；

2)为乘积r申请内存空间，长度为rn，并由指针rp指向，其中rn＝un+vn；

3)申请三个长度为rn的64位无符号长整型数组，分别由指针lpl，hpl，cl指向；

4)循环计算乘积的中间结果，缓存在步骤3)申请的三个临时数组中，其中中间结果的 每一位根据公式$\left(\mathrm{lpl}\right[\mathrm{ir}],\mathrm{cl}[\mathrm{ir}],\mathrm{hpl}[\mathrm{ir}\left]\right)=\underset{\mathrm{iu}+\mathrm{iv}=\mathrm{ir}}{\Sigma }u{p}_{\mathrm{iu}}*v{p}_{\mathrm{iv}}(0\le \mathrm{ir}<\mathrm{rn})$计算，lpl[ir]，cl[ir]和hpl[ir] 分别存储中间结果第ir位的低64位，中间64位以及高64位，且依据如下方式 循环计算：(rax,rdx)＝up_iu*vp_iv，(lpl[ir],c1)＝rax+lpl[ir]，(cl[ir],c2)＝cl[ir]+hpl[ir]+c1， hpl[ir]＝hpl[ir]+c2，此循环体采用内嵌汇编的方式实现，rax和rdx分别存储up_iu*vp_iv的低 64位和高64位，其中up_iu和vp_iv分别表示乘数u的第iu位和被乘数v第iv位，c1和c2表示进 位；

5)由步骤4)所述$\left(\mathrm{lpl}\right[\mathrm{ir}],\mathrm{cl}[\mathrm{ir}],\mathrm{hpl}[\mathrm{ir}\left]\right)=\underset{\mathrm{iu}+\mathrm{iv}=\mathrm{ir}}{\Sigma }u{p}_{\mathrm{iu}}*v{p}_{\mathrm{iv}}(0\le \mathrm{ir}<\mathrm{rn})$的循环计算过程采用 OpenMP的for任务分担策略进行并行化，开启多个计算线程并行计算，每个线程每次计算 由步骤4)所述的中间结果一个数位的计算流程，线程数可根据当前计算平台可用计算资源 决定，且给每个线程分配计算任务即中间结果数位的计算时采用动态调度策略，其中任务分 割粒度根据开启的线程数目以及计算机资源通过实验获得。如图2、3所示，本实验中将中 间结果的数位分配给8线程并行求取，数据分割粒度为rn/64，即将数位计算任务共分为 组，8个线程并行执行，每个线程每次领取rn/64个数位计算任务，率先执行 完任务的线程接着领取下组任务，无需等待；

6)对步骤4)和步骤5)所得的lpl[rn]，cl[rn]以及hpl[rn]三个数组中存储的中间结果 通过进位操作将最终结果存入rp[rn]数组中，设rlpl，rcl和rhpl存储进位数，则rp[rn]的第 ir位元素rp[ir](0≤ir<rn)以及进位依据如下方式计算：(rlpl,c1)＝rlpl+lpl[ir]， (rcl,c2)＝rcl+cl[ir]+c1，rhpl＝rhpl+hpl[ir]+c2，rp[ir]＝rlpl，rlpl＝rcl，rcl＝rhpl， rhpl＝0，其中(rlpl,c1)＝rlpl+lpl[ir]，(rcl,c2)＝rcl+cl[ir]+c1，rhpl＝rhpl+hpl[ir]+c2此 三个计算步骤采用内嵌汇编的方式实现，c1和c2表示进位；

7)计算乘积r的数位长度，如果rp_un+vn-1！＝0，则其长度为un+vn；否则，其长度为 un+vn-1。

测试结果如说明书附图中图4所示，测试平台为8核Intel Xeon 55系列处理器，CPU主 频为2.67GHz，以C语言为开发语言，Intel C++Compiler 11.1为编译器，并基于大整数运算 开源库GMP5.1.3版本进行算法的正确性验证和效率测试，测试时并行Comba算法为8线程， 其中数据规模为64位长整型双字为单位，程序运行时间单位为微秒。

本发明测试选取了10组大整数进行实验，每组数据为几百到上万数据单位的两个整数， 为了使得实验具有一般性，所有数据均是由程序随机产生。

实验中使用本发明中提出的并行执行程序与串行Comba算法在Intel的X86平台上进行 实验和性能对比。

从图4中可以看出，本发明选取的10组大整数进行实验得到的加速性能图，不同数据 规模下并行Comba算法与串行Comba算法的性能加速有所不同，其并行Comba算法与串行 Comba算法的加速比平均可达到4.30，性能提升了很多。

提供以上实施例仅仅是为了描述本发明的目的，而并非要限制本发明的范围。本发明的 范围由所附权利要求限定。不脱离本发明的精神和原理而做出的各种等同替换和修改，均应 涵盖在本发明的范围之内。